In [25]:
import pymc as pm
import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline
import pandas as pd
In [26]:
exp=np.genfromtxt("expdata.txt", delimiter='\t')
amp0=np.genfromtxt("amp0.txt", delimiter='\t')
iptg=exp[:,0][:,None]
tet0=exp[:,1][:,None]
MIC_exp=exp[:,2][:,None]
IND_exp=exp[:,3][:,None]
amp0=np.sort(amp0)[:,None]
In [27]:
constants_prob=np.genfromtxt("constants.txt", delimiter='\t')
In [28]:
def bandpass(constants_prob,iptg,tet0,amp0):
bg_bla=constants_prob[0]*5000 #20/0.1 #uM/min combo of beta and gamma
laci=constants_prob[1]*50 #0.01 #uM combo of laci/kd_laci
#kd_laci=constants_prob[2]*5e-3 #7.8e-4 #uM
n_laci=constants_prob[2]*5 #2.5
kd_iptg=constants_prob[3]*50 #25 #uM
n_iptg=constants_prob[4]*5 #2
#kcat=constants_prob[6]*5e6 #2.82e5 #1/m
km=constants_prob[5]*500 #100 #uM
ktr=constants_prob[6]*1e-2 #1e3 Combo of ktr and kcat
MIC_int=constants_prob[7]*5 #amp[0,5]
bg_tetc=constants_prob[8]*1000 #10/0.1 #uM/min combo of beta and gamma
ampr=constants_prob[9]*100 #0.1/1e-2 #uM combo of ampr/kd_ampr
n_ampr=constants_prob[10]*5 #1
kd_amp=constants_prob[11]*0.5 #0.5 #uM
n_amp=constants_prob[12]*5 #1
k_tet=constants_prob[13]*10 #10
MIC_tet_int=constants_prob[14]*5 #2
bla=np.zeros(len(iptg))
for i in range(len(iptg)):
D_iptg=1+(iptg[i]/kd_iptg)**n_iptg
D_laci=1+(laci/(D_iptg))**n_laci
bla[i]=bg_bla/D_laci
amp=np.zeros((len(bla),len(amp0)))
for i in range(len(bla)):
for j in range(len(amp0)):
coeff=[-ktr,(ktr*amp0[j]-ktr*km-bla[i]),ktr*amp0[j]*km]
sol=np.roots(coeff)
if sol[0]>=0 and sol[0]<=amp0[j]:
amp[i,j]=sol[0]
else:
amp[i,j]=sol[1]
MIC_mod=np.zeros(len(iptg))
for i in range(len(iptg)):
for j in range(len(amp0)):
if amp[i,j]>=MIC_int:
break
MIC_mod[i]=amp0[j]
tetc=np.zeros((len(iptg),len(amp0)))
tet=np.zeros((len(iptg),len(amp0)))
for i in range(len(iptg)):
for j in range(len(amp0)):
D_amp=1+(amp[i,j]/kd_amp)**n_amp
D_ampr=1+(ampr/(D_amp))**n_ampr
tetc[i,j]=bg_tetc/D_ampr
tet[i,j]=tet0[i]*k_tet/(k_tet+tetc[i,j])
#tet=tet0*k_tet_dif/(k_tet_dif+k_tet_act*tetc)
MIC_tet_mod=np.zeros(len(iptg))[:,None]
for i in range(len(iptg)):
for j in range(len(amp0)):
if tet[i,j]<=MIC_tet_int:
break
MIC_tet_mod[i]=amp0[j]
modeloutput=np.zeros((len(iptg),2))
for i in range(len(iptg)):
modeloutput[i,0]=MIC_mod[i]
modeloutput[i,1]=MIC_tet_mod[i]
return modeloutput
In [46]:
high_mod=np.zeros((constants_prob.shape[0],12))
low_mod=np.zeros((constants_prob.shape[0],12))
for i in range(constants_prob.shape[0]):
modeloutput=bandpass(constants_prob[i],iptg,tet0,amp0)
high_mod[i,:]=modeloutput[:,0]
low_mod[i,:]=modeloutput[:,1]
out=np.hstack((iptg,tet0,MIC_exp,high_mod.mean(axis=0)[:,None],IND_exp,low_mod.mean(axis=0)[:,None]))
pd.DataFrame(out,columns=['IPTG','Tet','Highpass(Exp)','Highpass(Model)','Lowpass(Exp)','Lowpass(Model)'])
Out[46]:
In [44]:
mean=high_mod.mean(axis=0)[:,None]
mean.shape
Out[44]:
In [ ]: