Upload New File

Ex_4_5.py

+#!/usr/bin/env python3
+# -*- coding: utf-8 -*-
+"""
+Created on Wed May  8 07:37:13 2024
+
+@author: delvigne
+"""
+
+'''
+Let's go back to exercice 3.3 with deterministic simulation of protein synthesis
+'''
+
+import numpy as np
+import matplotlib.pyplot as plt
+from scipy.integrate import odeint
+import random
+
+
+def SimpleGeneExpression(X,t):
+        alpha  = 0.09  #Removal rate in h-1
+        beta = 1.4   #Production rate in h-1
+        dXdt = beta - alpha*X
+        return dXdt
+
+t0   = 0
+tmax = 50
+dt   = 0.1
+t = np.arange(t0,tmax,dt)
+
+X0   = 0 #Initial protein concentration (g/L)
+X  = odeint(SimpleGeneExpression,X0,t)
+
+plt.plot(t,X)
+plt.xlabel('Time (h)')
+plt.ylabel('Protein X')
+plt.show()
+
+'''
+Now, let's turn it into a stochastic simulation
+'''
+X = [0]
+t = [0]
+tend = tmax
+
+alpha  = 0.09  #Removal rate in h-1
+beta = 1.4   #Production rate in h-1
+
+
+while t[-1] < tend:
+
+    props = [beta, alpha*X[-1]]
+
+
+    prop_sum = sum(props)
+
+
+    tau = np.random.exponential(scale=1/prop_sum)
+
+    t.append(t[-1] + tau)
+
+
+    rand = random.uniform(0,1)
+
+
+    # X production event
+    if rand * prop_sum <= props[0]:
+            X.append(X[-1] + 1)
+            
+
+    # X decay event
+    elif rand * prop_sum > props[0]:
+            X.append(X[-1] - 1)
+            
+
+plt.plot(t,X)
+plt.xlabel('Time (h)')
+plt.ylabel('Protein X')
+plt.show()
+