scripts.step_7_design_selection_and_map

step_7_design_selection_and_map.py

This script analyzes the simulated design space of genetic constructs to identify top-performing designs based on dynamic range and response time. It ranks designs using a composite score that favors high dynamic range and low response time (t50). The top candidates are visualized in a performance map alongside real experimental constructs for context. The final selected designs are saved for further experimental validation.


compute_metrics

compute_metrics(ts)

Compute performance metrics from the time series data.

Parameters:
  • ts (DataFrame) –

    Time series data with columns 'time' and 'Y'.

Returns: dict: Computed metrics: dynamic_range, t50, overshoot.

Source code in scripts/step_7_design_selection_and_map.py
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
def compute_metrics(ts):
    """
    Compute performance metrics from the time series data.

    Args:
        ts (pd.DataFrame): Time series data with columns 'time' and 'Y'.
    Returns:
        dict: Computed metrics: dynamic_range, t50, overshoot.
    """
    if ts is None or len(ts) < 5:
        return {k: np.nan for k in ["dynamic_range", "t50", "overshoot"]}

    y = ts["Y"].to_numpy()
    t = ts["time"].to_numpy()

    # 1. Dynamic Range
    dyn_range = y.max() - y.min()

    # 2. t50 (Time to 50% change)
    y0 = y[0]
    y_ss = y[-1]
    total_change = y_ss - y0
    if abs(total_change) < 1e-6:
        t50 = np.nan
    else:
        target = y0 + 0.5 * total_change
        idx = np.where(np.sign(y - target) != np.sign(y[0] - target))[0]
        if len(idx) > 0:
            i = idx[0]
            t1, t2 = t[i - 1], t[i]
            y1, y2 = y[i - 1], y[i]
            if y2 != y1:
                t50 = t1 + (target - y1) * (t2 - t1) / (y2 - y1)
            else:
                t50 = t1
        else:
            t50 = np.nan

    # 3. Overshoot
    tail_mean = y[int(0.8 * len(y)):].mean()
    overshoot = y.max() - tail_mean

    return {
        "dynamic_range": dyn_range,
        "t50": t50,
        "overshoot": overshoot
    }

ode_rhs

ode_rhs(t, y, t_up, K, n, alpha)

ODE right-hand side for the genetic construct model.

Args t (float): Current time. y (list): Current state [R, Y]. t_up (float): Half-time for upregulation. K (float): Hill constant. n (float): Hill coefficient. alpha (float): Degradation rate of Y.

Returns list: Derivatives [dR/dt, dY/dt].

Source code in scripts/step_7_design_selection_and_map.py
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
def ode_rhs(t, y, t_up, K, n, alpha):
    """
    ODE right-hand side for the genetic construct model.

    Args
        t (float): Current time.
        y (list): Current state [R, Y].
        t_up (float): Half-time for upregulation.
        K (float): Hill constant.
        n (float): Hill coefficient.
        alpha (float): Degradation rate of Y.

    Returns
        list: Derivatives [dR/dt, dY/dt].
    """
    R, Y = y
    t_up = max(float(t_up), 1e-6)
    beta = math.log(2.0) / t_up

    R_pos = max(R, 0.0)
    Kn = float(K) ** n
    Rn = R_pos ** n

    dR = beta - R_pos * (math.log(2.0) / t_up)
    dY = (Kn / (Kn + Rn)) - alpha * Y
    return [dR, dY]

simulate_construct

simulate_construct(t_up, K, n, alpha, t_end=72.0, dt=0.05)

Simulate the genetic construct dynamics over time.

Parameters:
  • t_up (float) –

    Half-time for upregulation.

  • K (float) –

    Hill constant.

  • n (float) –

    Hill coefficient.

  • alpha (float) –

    Degradation rate of Y.

  • t_end (float, default: 72.0 ) –

    End time for simulation.

  • dt (float, default: 0.05 ) –

    Time step for evaluation.

Returns: pd.DataFrame: Time series of Y over time.

Source code in scripts/step_7_design_selection_and_map.py
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
def simulate_construct(t_up, K, n, alpha, t_end=72.0, dt=0.05):
    """
    Simulate the genetic construct dynamics over time.

    Args:
        t_up (float): Half-time for upregulation.
        K (float): Hill constant.
        n (float): Hill coefficient.
        alpha (float): Degradation rate of Y.
        t_end (float): End time for simulation.
        dt (float): Time step for evaluation.
    Returns:
        pd.DataFrame: Time series of Y over time.
    """
    t_eval = np.arange(0.0, t_end + dt / 2, dt)
    y0 = [0.0, 1.0 / max(float(alpha), 1e-12)]

    sol = solve_ivp(
        fun=lambda t, y: ode_rhs(t, y, t_up, K, n, alpha),
        t_span=(0, t_end),
        y0=y0,
        t_eval=t_eval,
        method="LSODA",
        rtol=1e-7,
        atol=1e-9
    )
    if not sol.success:
        return None

    Y = sol.y[1] * alpha
    return pd.DataFrame({"time": t_eval, "Y": Y})