In [1]:
%pylab inline
import tellurium as te
from SloppyCell.ReactionNetworks import Utility
from scipy.stats import linregress

Populating the interactive namespace from numpy and matplotlib
In [2]:
#For reading excel files used in the modelling
!pip install xlwt

Requirement already satisfied: xlwt in /usr/local/lib/python2.7/dist-packages
In [3]:
import os
from pandas import *
import matplotlib.pyplot
from scipy.integrate import odeint
from lmfit import minimize, Parameter, Parameters, report_fit
from scipy.interpolate import interp1d, CubicSpline

##This allows is for writing in excel the results

from tempfile import TemporaryFile
from xlwt import Workbook

book = Workbook()

##########################

%run simple_model.py
%run parameters.py
In [4]:
def residuals(params, data, inter):
    y = array([params['y']])
    prediction=odeint(simple_model,y,range(0,24*7+2,1),args=(u1,params['degradation'], inter))
    p = {}
    for idx, i in enumerate(range(0,24*7+2,1)):
        p[i] = prediction[idx][0]

    res = []
    B_num = 0
    B_den = 0
    for k  in data.keys():
        B_num+=(p[k]*data[k])
        B_den+= p[k]**2
    B = B_num/B_den



    for k  in data.keys():
        res.append((B*p[k]-data[k])**2)
    
    return res

def sf(params, data, inter):
    y = array([params['y']])
    prediction=odeint(simple_model,y,range(0,24*7+2,1),args=(u1,params['degradation'], inter))
    p = {}
    for idx, i in enumerate(range(0,24*7+2,1)):
        p[i] = prediction[idx][0]

    res = []
    B_num = 0
    B_den = 0
    for k  in data.keys():
        B_num+=(p[k]*data[k])
        B_den+= p[k]**2
    B = B_num/B_den



    for k  in data.keys():
        res.append((B*p[k]-data[k])**2)
    
    return B

The RNA data is loaded first. I use the standard that I developed for my PhD. First we have the RNA data loaded into sloppy cell

In [5]:
book = ExcelFile('PhotoPeriodTiMetEdited.xlsx')
book.sheet_names

data = {}
for network in book.sheet_names:
    data[network] =book.parse(network)

genes = {'cL_m':'lhy','cC_m':'cca1','cP9_m':'prr9','cP7_m':'prr7','cP5_m':'prr5',
         'cT_m':'toc1','cLUX_m':'lux','cG_m':'gi','cE3_m':'elf3','cE4_m':'elf4'}

genotypes = {}

##From here we only use a subset of the data we are interested in constant light

#genotypes['pp6-18'] = []
#genotypes['pp8-16'] = []
#genotypes['pp12-12'] = []
#genotypes['pp18-6'] = []
genotypes['col_0_LL'] = []
#genotypes['col_0_DD'] = []
genotypes['prr79_LL'] = ['cP9_m','cP7_m']
#genotypes['prr79_DD'] = ['cP9_m','cP7_m']
genotypes['lhycca1_LL'] = ['cL_m']
#genotypes['lhycca1_DD'] = ['cL_m']
genotypes['toc1_LL'] = ['cT_m']
#genotypes['toc1_DD'] = ['cT_m']
genotypes['gi_LL'] = ['cG_m']
#genotypes['gi_DD'] = ['cG_m']


nuclei_number_gFW = 25.0e6
uncertainty = 0.6
displacement=24*10


## This function will be introduced in SloppyCell Utility module cos it is a very in order to make script more compact

def data_formater(data,genotypes,genes ,nuclei_number_gFW, displacement):
    '''
    Formating of data for sloppycell starting from an excel file
    This formating assumes that it is being fed with a dictionary of panda Data frames

    '''

    network_data={}

    for network in data.keys():
        if network not in genotypes.keys():
            continue

        gene_dict = {}
        for gene in genes.keys():
            try:
                if gene in genotypes[network]:
                    continue
            except:
                pass
            time_dict={}
            for zt in range(0,len(data[network][genes[gene]]),2):

                if not (np.isnan(data[network][genes[gene]][zt]) or np.isnan(data[network][genes[gene]][zt+1])):
                    points  = np.log(np.array([data[network][genes[gene]][zt],
                                               data[network][genes[gene]][zt+1]])/nuclei_number_gFW)
                    data_mean=np.mean(points)
                    data_var=np.std(points)                    
                else:
                    if np.isnan(data[network][genes[gene]][zt+1]) and not np.isnan(data[network][genes[gene]][zt]):
                        points  = np.log(np.array([data[network][genes[gene]][zt]])/nuclei_number_gFW)
                        data_mean=np.mean(points)
                        data_var = 0.6
                    elif np.isnan(data[network][genes[gene]][zt]) and not np.isnan(data[network][genes[gene]][zt+1]):
                        points  = np.log(np.array([data[network][genes[gene]][zt+1]])/nuclei_number_gFW)
                        data_mean=np.mean(points)
                        data_var = 0.6



                time_dict[zt+displacement] = (data_mean,data_var)
            gene_dict['log_'+gene]=time_dict
        network_data[network]=gene_dict
    return network_data

network_data = data_formater(data, genotypes,genes,nuclei_number_gFW, displacement)
In [6]:
exp(network_data['col_0_LL']['log_cL_m'][240][0])/exp(network_data['col_0_LL']['log_cC_m'][240][0])
Out[6]:
3.3138024885658273
In [7]:
figure(dpi=600)
def interpol_data(network_data, network, var):
    temp_time = array(sorted(network_data[network][var].keys()))
    temp = []
    for i in temp_time:
        temp.append(network_data[network][var][i][0])
    temp = array(temp)
    return temp_time-24*10-24, temp

interpolations = {}
for i in network_data['col_0_LL'].keys():
    temp_time, temp = interpol_data(network_data, 'col_0_LL', i)
    interpolations[i] = CubicSpline(temp_time,temp)

selected_genes = ['cL_m', 'cC_m']
colours = ['red','blue','orange','green']
k = 0
for j in selected_genes:
    plot(linspace(0,24*3,2000)-24,(interpolations['log_'+j](linspace(0,24*3,2000)-24)), lw=2,color=colours[k], label=genes[j].upper()+' mRNA')
    for i in array(sorted(network_data['col_0_LL']['log_'+j].keys())):
        errorbar(i-24*10-24, network_data['col_0_LL']['log_'+j][i][0], yerr=network_data['col_0_LL']['log_'+j][i][1], lw=2,fmt='s', color=colours[k])
    k+=1
ylim(-4,9)
xlim(-24,48)
xticks(range(0,48,12), fontsize=16)
yticks(fontsize=16)
ylabel('ln(#transcripts/cell)', fontsize=22)
xlabel('Time in constant light (h)', fontsize=22)
legend(loc='lower right', fontsize=16)
axvspan(-12,0, color='gray', alpha=0.7)
for i in range(0,3,1):
    axvspan(12+24*i,24*(i+1), color='gray', alpha=0.1)

tight_layout()
#savefig('cubic_interpolation_example.svg', format='svg', dpi=600, transparent=True)

Interpolation of protein CCA1 in entrainement conditions

In [8]:
figure(dpi=600)
temp_time, temp = interpol_data(network_data, 'col_0_LL', 'log_cL_m')
b = tile(temp[0:12],5)
b = concatenate((b,temp[12:]))
lhyinter = CubicSpline(range(0,len(b)*2,2),b)
plot(range(0,len(b)*2,2),exp(b))
plot(linspace(2,24*7,200),exp(lhyinter(linspace(2,24*7,200))))
xlim(24*3,24*7)
xlabel('Time in constant light (h)', fontsize=22)
axvspan(-12,0, color='gray', alpha=0.7)
for i in range(0,3,1):
    axvspan(12+24*i,24*(i+1), color='gray', alpha=0.1)
In [9]:
temp_time, temp = interpol_data(network_data, 'col_0_LL', 'log_cC_m')
b = tile(temp[0:12],5)
b = concatenate((b,temp[12:]))
cca1inter = CubicSpline(range(0,len(b)*2,2),b)
plot(range(0,len(b)*2,2),exp(b))
plot(linspace(2,24*7,200),exp(cca1inter(linspace(2,24*7,200))))
Out[9]:
[<matplotlib.lines.Line2D at 0x4081d9f910>]

Simple model that has one differential equation. We do not fitt the model to the data as there are degradation rates for LHY. The degradation rate comes from Hansen Louise un published

In [10]:
from scipy.integrate import odeint
def simple_model(y,t,u1,d1,mdata1):
    mRNA = mdata1(t)
    Prot0 = y[0]
    dydt=range(1)
    dydt[0]= u1*exp(mRNA)-d1*Prot0
    return dydt
In [11]:
#in aminoacids 
amin_sec = 3.
rib_tranc = 10.
CCA1_length = 609.
CCA1_trans_rate = 1/(CCA1_length/amin_sec/10)*3600
print CCA1_trans_rate
amin_sec = 3.
rib_tranc = 10.
LHY_length = 646.
LHY_trans_rate = 1/(LHY_length/amin_sec/10)*3600
print LHY_trans_rate

177.339901478
167.182662539

Model behaviour in 12L:12D using the last LD cycle of Wang figure 1

In [16]:
figure(dpi=600)
deg_cca1 = read_excel('CCA1_degradation_louise.xlsx', sheet_name='Sheet1')
CCA1_reg = linregress(deg_cca1.iloc[:,0]/60.,log(mean(deg_cca1.iloc[:,1:4], axis=1)/mean(deg_cca1.iloc[:,1:4], axis=1)[0]))
plot(deg_cca1.iloc[:,0]/60.,log(mean(deg_cca1.iloc[:,1:4], axis=1)/mean(deg_cca1.iloc[:,1:4], axis=1)[0]),"-o", color="black",label="data")
plot(deg_cca1.iloc[:,0]/60., CCA1_reg.slope*deg_cca1.iloc[:,0]/60.+CCA1_reg.intercept, color="black", label="reg.")
xlabel('Time in hours', fontsize=12)
ylabel('ln(CCA1 signal/(CCA1 t0))', fontsize=12)
legend(loc="upper right")
tight_layout()
savefig("CCA1_deacy_linear_reg_hansen.svg", format="svg", dpi=600)
savefig("CCA1_deacy_linear_reg_hansen.pdf", format="pdf", dpi=600)
In [17]:
#the second sheet has the data for LHY from the lhy-1 genotype in which the protein is overexpressed with
#transposone this is from Hae-Ryong Song and Carre 2005
deg_cca1 = read_excel('CCA1_degradation_louise.xlsx', sheet_name='Sheet2')
lhy_reg = linregress(deg_cca1.time/60,deg_cca1.deg)
lhy_reg.slope
Out[17]:
-1.7688
In [38]:
pCCA1LD = read_excel('ProteinTimeseriesDB2.xlsx', sheet_name='pCCA1(wang1998LD)')
plot(pCCA1LD.time,pCCA1LD.pCCA1LD/max(pCCA1LD.pCCA1LD),'s-',color="blue")
pLHY = read_excel('ProteinTimeseriesDB2.xlsx',sheet_name='pLHY(karre)')
plot(pLHY.time,pLHY.pLHYLD/100,"s-",color="red")
xticks(range(0,24*6,12))
pCCA1 = read_excel('ProteinTimeseriesDB2.xlsx', sheet_name='pCCA1(Wang1998)')
plot(pCCA1.time+72,pCCA1.pCCA1,'o--',color="blue")
Out[38]:
[<matplotlib.lines.Line2D at 0x40aee55890>]
In [19]:
res = {}
pro_data = {}
In [20]:
#this has units of hour, you can see that lhy is double the stability of CCA1. 
print log(2)/lhy_reg.slope*(-1)
print log(2)/CCA1_reg.slope*(-1)

0.3918742540479112
0.20434558391843624
In [48]:
figure(dpi=600)
y = zeros(1)
t = range(0,168,1)
t2 = numpy.array(range(1,168,1))
res['LHY'] = odeint(simple_model,y,t,args=(LHY_trans_rate,lhy_reg.slope*(-1),lhyinter))
#plot(pCCA1.time,log(pCCA1.pCCA1*10000), 'bo-', label='Knowles2008')

xticks(range(0,24*4,12))
plot(t2-24*5,res['LHY'][1:,0], lw=3, color= 'red', label='pLHY')
yticks(fontsize=16)
xticks(fontsize=16)
ylabel('#monomers/cell', fontsize=22)
xlabel('Time in h', fontsize=22)
xticks(fontsize=16)
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(24+12,24*2, color='gray',alpha=0.3)
axvspan(24*2+12,24*3, color='gray',alpha=0.3)
axvspan(12,24, color='gray',alpha=0.3)
y = zeros(1)
res['CCA1'] = odeint(simple_model,y,t,args=(CCA1_trans_rate,CCA1_reg.slope*(-1),cca1inter))
plot(t2-24*5,(res['CCA1'][1:,0]), lw=3, color= 'blue', label='pCCA1')
#plot(pCCA1.time,pCCA1.pCCA1*10000,"-o", color="blue")
xlim(-24,24*2-4+24)
yticks(fontsize=16)
xticks(fontsize=16)
ylim(1,1e6)
ylabel('#monomers/cell', fontsize=22)
xlabel('Time in constant light (h)', fontsize=22)
xticks(fontsize=16)
plot(t2-24*5,(res['LHY'][1:,0]+res['CCA1'][1:,0]), lw=3, color= 'orange', label='pLHY+pCCA1')

pro_data['cL'] = {}
for idx, i in enumerate(res['CCA1']):
    if idx >= 24*4:
        pro_data['cL'][t[idx]] =  (res['CCA1'][idx][0]+res['LHY'][idx][0],0.1)


pCCA1LD = read_excel('ProteinTimeseriesDB2.xlsx', sheet_name='pCCA1(wang1998LD)')
plot(pCCA1LD.time-72,pCCA1LD.pCCA1LD/max(pCCA1LD.pCCA1LD)*10000,'s-',color="blue", label="CCA1 LD")
pLHY = read_excel('ProteinTimeseriesDB2.xlsx',sheet_name='pLHY(karre)')
plot(pLHY.time-72+24,pLHY.pLHYLD/100*100000,"s-",color="red", label="LHY LD")
xticks(range(0,24*3,12))
pCCA1 = read_excel('ProteinTimeseriesDB2.xlsx', sheet_name='pCCA1(Wang1998)')
plot(pCCA1.time,pCCA1.pCCA1*10000,'o--',color="blue", label="CCA1 LL")
        
        
yscale('log')
legend(loc='lower right', fontsize=16)
tight_layout()
#savefig('lhy_cca1_simple_model.svg', format='svg', dpi=600, transparent=True)
savefig('lhy_cca1_simple_model_plus_data.svg', format='svg', dpi=600, transparent=True)
savefig('lhy_cca1_simple_model_plus_data.pdf', format='pdf', dpi=600, transparent=True)
In [49]:
def simple_model_ligh_deg(y,t,u1,d1,d2,period,dusk,tw,mdata1):
    if t < 24*5:
        l=0.5*((1+np.tanh((t-period*np.floor(t/period))/tw))-(1+np.tanh((t-period*np.floor(t/period)-dusk)/tw))+(1+np.tanh((t-period*np.floor(t/period)-period)/tw)))
    else :
        l=1.0
    mRNA = mdata1(t)
    Prot0 = y[0]
    dydt=range(1)
    dydt[0]= u1*exp(mRNA)-(d1+d2*(1-l))*Prot0
    return dydt
In [50]:
#This comes from Eva Farre 2007
figure(dpi=600)
prr7decay = read_excel('data/prr7prr9_decay.xlsx',sheet_name='PRR7decaycurve')
     
PRR7_reg_light = linregress(prr7decay.timeH,log(prr7decay.ZT12L))
plot(prr7decay.timeH,log(prr7decay.ZT12L), "o-", color="blue", label="light")
plot(prr7decay.timeH, PRR7_reg_light.slope*prr7decay.timeH+PRR7_reg_light.intercept,color="blue", label="reg.")

PRR7_reg_dark = linregress(prr7decay.timeH,log(prr7decay.ZT12D))
plot(prr7decay.timeH,log(prr7decay.ZT12D), "o-", color="black", label="dark")
plot(prr7decay.timeH, PRR7_reg_dark.slope*prr7decay.timeH+PRR7_reg_dark.intercept, color="black", label="reg.")
xlabel("Time in hours", fontsize=12)
ylabel("ln(normalised decay)", fontsize=12)
legend(loc="upper right")
tight_layout()
savefig("PRR7_decay_constant_regression_farre.svg", format="svg", dpi=600)
savefig("PRR7_decay_constant_regression_farre.pdf", format="pdf", dpi=600)
In [51]:
temp_time, temp = interpol_data(network_data, 'col_0_LL', 'log_cP7_m')
b = tile(temp[0:12],5)
b = concatenate((b,temp[12:]))
p7inter = CubicSpline(range(0,len(b)*2,2),b)
plot(range(0,len(b)*2,2),exp(b))
plot(linspace(2,24*7,200),exp(p7inter(linspace(2,24*7,200))))
Out[51]:
[<matplotlib.lines.Line2D at 0x40d0a38c90>]
In [52]:
figure(dpi=600)
amin_sec = 3.
rib_tranc = 10.
PRR7_length = 727.
PRR7_trans_rate = 1/(PRR7_length/amin_sec/10)*3600
print PRR7_trans_rate
y = zeros(1)
period = 24.
dusk = 12.
tw = 1.
daypat = 'LD'
t = numpy.array(range(0,24*7,1))
res['PRR7'] = odeint(simple_model_ligh_deg,y,t,
              args=(PRR7_trans_rate,PRR7_reg_light.slope*(-1),PRR7_reg_dark.slope*(-1)-PRR7_reg_light.slope*(-1),
                    period,dusk,tw,
                    p7inter))
pro_data['cP7'] = {}
for idx, i in enumerate(res['PRR7']):
    if idx >= 24*4:
        pro_data['cP7'][t[idx]] =  (res['PRR7'][idx][0],0.1)

        
plot(t-24*5,res['PRR7'][:,0], lw=3, color= 'black', label='LHY')
ylabel('#monomers/cell', fontsize=22)
xlabel('Time in h', fontsize=22)
xticks(fontsize=16)
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(12,24, color='gray',alpha=0.1)
axvspan(24+12,24*3, color='gray',alpha=0.1)
xticks(range(0,24*3,12),fontsize=16)
ylabel('#monomers/cell', fontsize=22)
xlabel('Time in contant light (h)', fontsize=22)
xlim(-24,24*2)
yticks(fontsize=16)

148.555708391
Out[52]:
(array([-10000.,      0.,  10000.,  20000.,  30000.,  40000.,  50000.,
         60000.,  70000.]), <a list of 9 Text yticklabel objects>)
In [28]:
print (pro_data['cP7'][96+9][0]-pro_data['cP7'][96+24][0])/2
print (pro_data['cP7'][96+9+24][0]-pro_data['cP7'][96+24*2+1][0])/2
print (pro_data['cP7'][96+9+24][0]-pro_data['cP7'][96+24*2+1][0])/(pro_data['cP7'][96+9][0]-pro_data['cP7'][96+24][0])

15411.987765454134
24538.23309287341
1.5921523859417859
In [55]:
#from ito shogo 2007 
figure(dpi=600)
prr7decay = read_excel('data/prr7prr9_decay.xlsx',sheet_name='PRR9decay')
     
PRR7_reg_light = linregress(prr7decay.timemin,log(prr7decay.morning/max(prr7decay.morning)))
plot(prr7decay.timemin,log(prr7decay.morning/max(prr7decay.morning)), "o-",color="blue")
plot(prr7decay.timemin, PRR7_reg_light.slope*prr7decay.timemin+PRR7_reg_light.intercept, color="blue")
ylabel("ln(normalised decay)", fontsize=12)
xlabel("Time in hours", fontsize=12)
tight_layout()
savefig("PRR9_decay_normalised_ito_shogo.svg", format="svg", dpi=600)
savefig("PRR9_decay_normalised_ito_shogo.pdf", format="svg", dpi=600)
In [56]:
temp_time, temp = interpol_data(network_data, 'col_0_LL', 'log_cP9_m')
b = tile(temp[0:12],5)
b = concatenate((b,temp[12:]))
p7inter = CubicSpline(range(0,len(b)*2,2),b)
plot(range(0,len(b)*2,2),exp(b))
plot(linspace(2,24*7,200),exp(p7inter(linspace(2,24*7,200))))
Out[56]:
[<matplotlib.lines.Line2D at 0x40cb969f10>]
In [57]:
figure(figsize=(10,7))
amin_sec = 3.
rib_tranc = 10.
PRR7_length = 468.
PRR7_trans_rate = 1/(PRR7_length/amin_sec/10)*3600
print PRR7_trans_rate
y = zeros(1)
period = 24.
dusk = 12.
tw = 1.
daypat = 'LD'
t = numpy.array(range(0,24*7,1))
res['PRR9'] = odeint(simple_model,y,t,
              args=(PRR7_trans_rate,PRR7_reg_light.slope*(-1),
                    p7inter))
plot(t-24*5,res['PRR9'][:,0], lw=3, color= 'black', label='LHY')
pro_data['cP9'] = {}
for idx, i in enumerate(res['PRR9']):
    if idx >= 24*4:
        pro_data['cP9'][t[idx]] =  (res['PRR9'][idx][0],0.1)

xticks(fontsize=16)
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(24+12,24*2, color='gray',alpha=0.1)
axvspan(24*2+12,24*3, color='gray',alpha=0.1)
axvspan(12,24, color='gray',alpha=0.1)
xticks(range(0,24*3,12),fontsize=16)
ylabel('#monomers/cell', fontsize=22)
xlabel('Time in contant light (h)', fontsize=22)
xlim(-24,24*2)
yticks(fontsize=16)
axvspan(12,24, color='gray',alpha=0.7)
axvline(2)
axvline(8)

230.769230769
Out[57]:
<matplotlib.lines.Line2D at 0x40a79c9a10>
In [58]:
#data from takoshi kiba 2007
temp_time, temp = interpol_data(network_data, 'col_0_LL', 'log_cP5_m')
b = tile(temp[0:12],5)
b = concatenate((b,temp[12:]))
p7inter = CubicSpline(range(0,len(b)*2,2),b)
plot(range(0,len(b)*2,2),exp(b))
plot(linspace(2,24*7,200),exp(p7inter(linspace(2,24*7,200))))

prr7decay = read_excel('data/prr7prr9_decay.xlsx',sheetname='PRR5decay')
     
PRR7_reg_wl = linregress(prr7decay.time,log(prr7decay.wl))
plot(prr7decay.time,log(prr7decay.wl))
plot(prr7decay.time, PRR7_reg_wl.slope*prr7decay.time+PRR7_reg_wl.intercept)

PRR7_reg_dark = linregress(prr7decay.time,log(prr7decay.dark))
plot(prr7decay.time,log(prr7decay.dark))
plot(prr7decay.time, PRR7_reg_dark.slope*prr7decay.time+PRR7_reg_dark.intercept)

figure(figsize=(10,7))
amin_sec = 3.
rib_tranc = 10.
PRR7_length = 558.
PRR7_trans_rate = 1/(PRR7_length/amin_sec/10)*3600
print PRR7_trans_rate
y = zeros(1)
period = 24.
dusk = 12.
tw = 1.
daypat = 'LD'
t = numpy.array(range(0,24*7,1))
res['PRR5'] = odeint(simple_model_ligh_deg,y,t,
              args=(PRR7_trans_rate,PRR7_reg_wl.slope*(-1),PRR7_reg_dark.slope*(-1)-PRR7_reg_wl.slope*(-1),
                    period,
                    dusk,
                    tw,
                    p7inter))

plot(t-24*5,res['PRR5'][:,0], lw=3, color= 'black', label='LHY')
pro_data['cP5'] = {}
for idx, i in enumerate(res['PRR5']):
    if idx >= 24*4:
        pro_data['cP5'][t[idx]] =  (res['PRR5'][idx][0],0.1)

xticks(fontsize=16)
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(24+12,24*2, color='gray',alpha=0.1)
axvspan(24*2+12,24*3, color='gray',alpha=0.1)
axvspan(12,24, color='gray',alpha=0.1)
xticks(range(0,24*3,12),fontsize=16)
ylabel('#monomers/cell', fontsize=22)
xlabel('Time in contant light (h)', fontsize=22)
xlim(-24,24*2)
yticks(fontsize=16)

/usr/local/lib/python2.7/dist-packages/pandas/util/_decorators.py:188: FutureWarning:

The `sheetname` keyword is deprecated, use `sheet_name` instead


193.548387097
Out[58]:
(array([-10000.,      0.,  10000.,  20000.,  30000.,  40000.,  50000.,
         60000.,  70000.]), <a list of 9 Text yticklabel objects>)
In [59]:
#TOC1 degradation comes from Paloma mas Science paper
temp_time, temp = interpol_data(network_data, 'col_0_LL', 'log_cT_m')
b = tile(temp[0:12],5)
b = concatenate((b,temp[12:]))
p7inter = CubicSpline(range(0,len(b)*2,2),b)
plot(range(0,len(b)*2,2),exp(b))
plot(linspace(2,24*7,200),exp(p7inter(linspace(2,24*7,200))))
figure(figsize=(10,7))
prr7decay = read_excel('data/prr7prr9_decay.xlsx',sheetname='TOC1decay')
     
PRR7_reg_wl = linregress(prr7decay.time,log(prr7decay.normdecay))
plot(prr7decay.time,log(prr7decay.normdecay))
plot(prr7decay.time, PRR7_reg_wl.slope*prr7decay.time+PRR7_reg_wl.intercept)

figure(figsize=(10,7))
amin_sec = 3.
rib_tranc = 10.
PRR7_length = 618.
PRR7_trans_rate = 1/(PRR7_length/amin_sec/10)*3600
print PRR7_trans_rate
y = zeros(1)
period = 24.
dusk = 12.
tw = 1.
daypat = 'LD'
t = numpy.array(range(0,24*7,1))
res['TOC1'] = odeint(simple_model,y,t,
              args=(PRR7_trans_rate,PRR7_reg_wl.slope*(-1),
                    p7inter))
plot(t-24*5,res['TOC1'][:,0], lw=3, color= 'black', label='TOC1')

pro_data['cT'] = {}
for idx, i in enumerate(res['TOC1']):
    if idx >= 24*4:
        pro_data['cT'][t[idx]] =  (res['TOC1'][idx][0],0.1)

xticks(fontsize=16)
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(24+12,24*2, color='gray',alpha=0.1)
axvspan(24*2+12,24*3, color='gray',alpha=0.1)
axvspan(12,24, color='gray',alpha=0.1)
xticks(range(0,24*3,12),fontsize=16)
ylabel('#monomers/cell', fontsize=22)
xlabel('Time in contant light (h)', fontsize=22)
xlim(-24,24*2)
yticks(fontsize=16)

174.757281553
Out[59]:
(array([-5000.,     0.,  5000., 10000., 15000., 20000., 25000., 30000.,
        35000.]), <a list of 9 Text yticklabel objects>)
In [61]:
figure(dpi=600)
t = numpy.array(range(0,24*7,1))
plot(t-24*5,res['PRR9'][:,0], lw=3, color= 'black', label='P9')
pProtein = read_excel('ProteinTimeseriesDB2.xlsx', sheetname='PRR9(Fujiware2008fig1F)')
plot(pProtein.time,pProtein.pPRR9*9000, 'ko-')
w_data={}
w_data['cP9'] ={}
for idx,i in enumerate(pProtein.time):
     w_data['cP9'][int(pProtein.iloc[idx,0])+24*5]  =  pProtein.iloc[idx,1]

plot(t-24*5,res['PRR7'][:,0], lw=3, color= 'blue', label='P7')
pProtein = read_excel('ProteinTimeseriesDB2.xlsx', sheetname='pPRR7(fujiware2008fig1F)')
plot(pProtein.time,pProtein.pPRR7*60000, 'bo-')
w_data['cP7'] ={}
for idx,i in enumerate(pProtein.time):
     w_data['cP7'][int(pProtein.iloc[idx,0])+24*5]  =  pProtein.iloc[idx,1]
xticks(fontsize=16)
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(24+12,24*2, color='gray',alpha=0.3)
axvspan(24*2+12,24*3, color='gray',alpha=0.3)
axvspan(12,24, color='gray',alpha=0.3)
xticks(range(0,24*3,12),fontsize=16)
ylabel('#monomers/cell', fontsize=22)
xlabel('Time in contant light (h)', fontsize=22)
xlim(-24,24*2)
yticks(fontsize=16)
ylim(1,1e6)
yscale('log')
legend(loc='lower left', fontsize=16)
tight_layout()
#savefig('PRR9_PRR7_protein_prediction.svg',format='svg',dpi=600,transparent=True)
#yscale('log')
In [63]:
figure(dpi=600)
plot(t-24*5,res['PRR5'][:,0], lw=3, color= 'orange', label='P5')
pProtein = read_excel('ProteinTimeseriesDB2.xlsx', sheetname='pPRR5(Fujiware2008fig1F)')
plot(pProtein.time,pProtein.pPRR5*60000, 'o-', color='orange')
w_data['cP5'] ={}
for idx,i in enumerate(pProtein.time):
     w_data['cP5'][int(pProtein.iloc[idx,0])+24*5]  =  pProtein.iloc[idx,1]

plot(t-24*5,res['TOC1'][:,0], lw=3, color= 'red', label='TOC1')
pProtein = read_excel('ProteinTimeseriesDB2.xlsx', sheetname='pTOC1(Fujiware2008fig1F)')
plot(pProtein.time,pProtein.pTOC1*30000, 'ro-')
w_data['cT'] ={}
for idx,i in enumerate(pProtein.time):
     w_data['cT'][int(pProtein.iloc[idx,0])+24*5]  =  pProtein.iloc[idx,1]
xticks(fontsize=16)
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(24+12,24*2, color='gray',alpha=0.3)
axvspan(24*2+12,24*3, color='gray',alpha=0.3)
axvspan(12,24, color='gray',alpha=0.3)
xticks(range(0,24*3,12),fontsize=16)
ylabel('#monomers/cell', fontsize=22)
xlabel('Time in contant light (h)', fontsize=22)
xlim(-24,24*2)
ylim(1,1e6)
yticks(fontsize=16)
legend(loc='lower right', fontsize=16)
yscale('log')
tight_layout()
#savefig('PRR5_TOC1_protein_prediction.svg',format='svg',dpi=600,transparent=True)

#savefig('figures/PRR51_simple_model.png', format='png', dpi=300)
In [64]:
def simple_model(y,t,u1,d1,mdata1):
    mRNA = mdata1(t)
    Prot0 = y[0]
    dydt=range(1)
    dydt[0]= u1*exp(mRNA)-d1*Prot0
    return dydt
In [65]:
def residuals(params, data, inter):
    y = array([params['y']])
    prediction=odeint(simple_model,y,range(0,24*7+2,1),args=(u1,params['degradation'], inter))
    p = {}
    for idx, i in enumerate(range(0,24*7+2,1)):
        p[i] = prediction[idx][0]

    res = []
    B_num = 0
    B_den = 0
    for k  in data.keys():
        B_num+=(p[k]*data[k])
        B_den+= p[k]**2
    B = B_num/B_den



    for k  in data.keys():
        res.append((B*p[k]-data[k])**2)
    
    return res

def sf(params, data, inter):
    y = array([params['y']])
    prediction=odeint(simple_model,y,range(0,24*7+2,1),args=(u1,params['degradation'], inter))
    p = {}
    for idx, i in enumerate(range(0,24*7+2,1)):
        p[i] = prediction[idx][0]

    res = []
    B_num = 0
    B_den = 0
    for k  in data.keys():
        B_num+=(p[k]*data[k])
        B_den+= p[k]**2
    B = B_num/B_den



    for k  in data.keys():
        res.append((B*p[k]-data[k])**2)
    
    return B
In [67]:
t = numpy.array(t)

pELF3 = read_excel('data/ELF3ELF4LUXProtein.xlsx', sheetname='pELF4')

w_data['cE4'] ={}
for idx,i in enumerate(pELF3.time):
     w_data['cE4'][int(pELF3.iloc[idx,0])+24*5]  =  pELF3.iloc[idx,1]

temp_time, temp = interpol_data(network_data, 'col_0_LL', 'log_cE4_m')
b = tile(temp[0:12],5)
b = concatenate((b,temp[12:]))
p7inter = CubicSpline(range(0,len(b)*2,2),b)
plot(range(0,len(b)*2,2),exp(b))
plot(linspace(2,24*7,200),exp(p7inter(linspace(2,24*7,200))))

figure(dpi=600)
xticks(fontsize=16)
xticks(range(0,24*3,12),fontsize=16)
param = Parameters()
param.add('y',value=1.0173e+05, min=0.0, max=1e6, vary=True)
param.add('degradation',value=0.3848250, min=0.0, max=3.0, vary=True)
out = minimize(residuals,param,args=(w_data['cE4'],p7inter))
report_fit(out.params)
res['ELF4']=odeint(simple_model,array([out.params['y'].value]),t,args=(u1,out.params['degradation'].value,p7inter))
B = sf(out.params, w_data['cE4'],p7inter)
plot(t-24*5,res['ELF4'], 'red', label='pELF4', lw=3)
plot(pELF3.time-24, pELF3.pELF4/B, 'ro-', label='ELF4')
xlabel('Time in constant light (h)', fontsize=22)
xlim(-24,24*2)
#ylim(0,1.4e5)
yticks(fontsize=16)
ylabel('#monomers/cell', fontsize=22)

pro_data['cE4'] = {}
for idx, i in enumerate(res['ELF4']):
    if idx >= 24*4:
        pro_data['cE4'][t[idx]] =  (res['ELF4'][idx][0],0.1)



pELF3 = read_excel('data/ELF3ELF4LUXProtein.xlsx', sheetname='pELF3')
w_data['cE3'] ={}
for idx,i in enumerate(pELF3.time):
     w_data['cE3'][int(pELF3.iloc[idx,0])+24*5]  =  pELF3.iloc[idx,1]

temp_time, temp = interpol_data(network_data, 'col_0_LL', 'log_cE3_m')
b = tile(temp[0:12],5)
b = concatenate((b,temp[12:]))
p7inter = CubicSpline(range(0,len(b)*2,2),b)
#plot(range(0,len(b)*2,2),exp(b))
#plot(linspace(2,24*7,200),exp(p7inter(linspace(2,24*7,200))))


xticks(fontsize=16)
xticks(range(0,24*3,12),fontsize=16)
param = Parameters()
param.add('y',value=1.0173e+05, min=0.0, max=1e6, vary=True)
param.add('degradation',value=0.30848250, min=0.0, max=3.0, vary=True)
out = minimize(residuals,param,args=(w_data['cE3'],p7inter))
report_fit(out.params)
res['ELF3']=odeint(simple_model,array([out.params['y'].value]),t,args=(u1,out.params['degradation'].value,p7inter))
B = sf(out.params, w_data['cE3'],p7inter)
plot(t-24*5,res['ELF3'], 'blue', label='pELF3', lw=3)
plot(pELF3.time-24, pELF3.pELF3/B, 'bo-', label='cE3')
xlabel('Time in constant light (h)', fontsize=22)
xlim(-24,24*2)
#ylim(0,30e3)
yticks(fontsize=16)
ylabel('#monomers/cell', fontsize=22)

pro_data['cE3'] = {}
for idx, i in enumerate(res['ELF3']):
    if idx >= 24*4:
        pro_data['cE3'][t[idx]] =  (res['ELF3'][idx][0],0.1)


pELF3 = read_excel('data/ELF3ELF4LUXProtein.xlsx', sheetname='pLUX')

w_data['cLUX'] ={}
for idx,i in enumerate(pELF3.time):
     w_data['cLUX'][int(pELF3.iloc[idx,0])+24*5]  =  pELF3.iloc[idx,1]

temp_time, temp = interpol_data(network_data, 'col_0_LL', 'log_cLUX_m')
b = tile(temp[0:12],5)
b = concatenate((b,temp[12:]))
p7inter = CubicSpline(range(0,len(b)*2,2),b)
#plot(range(0,len(b)*2,2),exp(b))
#plot(linspace(2,24*7,200),exp(p7inter(linspace(2,24*7,200))))

xticks(fontsize=16)
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(24+12,24*2, color='gray',alpha=0.3)
axvspan(12,24, color='gray',alpha=0.3)
xticks(range(0,24*3,12),fontsize=16)
param = Parameters()
param.add('y',value=1.0173e+05, min=0.0, max=1e6, vary=True)
param.add('degradation',value=0.03848250, min=0.0, max=3.0, vary=True)
out = minimize(residuals,param,args=(w_data['cLUX'],p7inter))
report_fit(out.params)
res['LUX']=odeint(simple_model,array([out.params['y'].value]),t,args=(u1,out.params['degradation'].value,p7inter))
B = sf(out.params, w_data['cLUX'],p7inter)

pro_data['cLUX'] = {}
for idx, i in enumerate(res['LUX']):
    if idx >= 24*4:
        print i
        pro_data['cLUX'][t[idx]] =  (res['LUX'][idx][0],0.1)


plot(t-24*5,res['LUX'], 'black', label='pLUX', lw=3)
plot(pELF3.time-24, pELF3.pLUX/B, 'ko-', label='LUX')
xlabel('Time in constant light (h)', fontsize=22)
xlim(-24,24*2)
ylim(1,1e6)
yticks(fontsize=16)
ylabel('#monomers/cell', fontsize=22)
legend(loc='upper right', fontsize=16)
yscale('log')
#savefig('EC_protein_prediction.svg',format='svg',dpi=600,transparent=True)

#savefig('figures/EC_fitting.png', format='png',dpi=300)

[[Variables]]
    y:            101729.934 +/- 18002.0879 (17.70%) (init = 101730)
    degradation:  0.38482477 +/- 0.04124358 (10.72%) (init = 0.384825)
[[Correlations]] (unreported correlations are < 0.100)
    C(y, degradation) = -0.929
[[Variables]]
    y:            107805.104 +/- 6407.37785 (5.94%) (init = 101730)
    degradation:  0.07779919 +/- 0.01510069 (19.41%) (init = 0.3084825)
[[Correlations]] (unreported correlations are < 0.100)
    C(y, degradation) = -0.922
[[Variables]]
    y:            103360.576 +/- 1330503.32 (1287.24%) (init = 101730)
    degradation:  0.04196282 +/- 0.01075398 (25.63%) (init = 0.0384825)
[[Correlations]] (unreported correlations are < 0.100)
    C(y, degradation) =  0.940
[70649.94604933]
[68403.92170012]
[66178.52141009]
[63977.54112552]
[61829.86885218]
[59975.76298822]
[59001.80089699]
[59936.15075344]
[63953.60521998]
[71752.10111323]
[82309.24772231]
[92289.83904187]
[98192.8756314]
[99230.63653349]
[97323.56537845]
[94321.21029063]
[90906.1685419]
[87314.08625754]
[83769.00770673]
[80380.76612246]
[77337.15153228]
[75153.01740683]
[73809.42111158]
[72273.03533385]
[70135.08181061]
[67742.24195775]
[65356.54687234]
[63097.75811601]
[61034.2622367]
[59179.61413707]
[57675.64761174]
[57535.76597543]
[62543.18526579]
[77558.65212044]
[96939.99927261]
[108419.59499932]
[111433.32975803]
[111452.84587023]
[111497.31868486]
[112969.39739504]
[115689.99863066]
[116841.44072113]
[115511.65084306]
[113409.02262869]
[111722.81337717]
[110464.07473148]
[108865.31122066]
[106364.46947618]
[103272.88519688]
[100351.24309622]
[98101.10694738]
[95881.41515551]
[93195.03433159]
[90735.27061503]
[89379.74504352]
[88749.8101671]
[88048.14539937]
[88928.40985501]
[95363.28395825]
[106528.23756302]
[114980.38845825]
[119369.24048392]
[121246.30176596]
[121469.80178563]
[120829.81396236]
[120406.33129398]
[120974.40300529]
[122020.29904494]
[122351.32696723]
[121455.01765552]
[119481.19470812]
[116835.91272031]

Here is were use the new models for rescaling¶

In [68]:
import tellurium as te
from SloppyCell.ReactionNetworks import Utility
In [73]:
%pylab inline
U2017 = te.loada('Models/U2019.4.txt')
U2017res = U2017.simulate(0,24*10,300)
plot(U2017res['time'],U2017res['[cL]']/max(U2017res['[cL]']))
U2017 = te.loada('Models/U2019.4.txt')
U2017.setValue('gl',5)
U2017.setValue('cL',U2017.getValue('cL')*U2017.getValue('gl'))
U2017res = U2017.simulate(0,24*10,300)
#plot(res['time'],res['[cL]'])

plot(U2017res['time'],U2017res['[cL]']/max(U2017res['[cL]']))

Populating the interactive namespace from numpy and matplotlib
Out[73]:
[<matplotlib.lines.Line2D at 0x40bc531290>]
In [74]:
U2017 = te.loada('Models/U2019.4.txt')
U2017res = U2017.simulate(0,24*10,300)
plot(U2017res['time'],U2017res['[cP7]']/max(U2017res['[cP7]']))
U2017 = te.loada('Models/U2019.4.txt')
U2017.setValue('gp7',2)
U2017.setValue('cP9',U2017.getValue('cP9')*U2017.getValue('gp7'))
U2017.setValue('cP7',U2017.getValue('cP7')*U2017.getValue('gp7'))
U2017.setValue('cP5',U2017.getValue('cP5')*U2017.getValue('gp7'))
U2017.setValue('cT',U2017.getValue('cT')*U2017.getValue('gp7'))
res = U2017.simulate(0,24*10,300)
plot(U2017res['time'],U2017res['[cP7]']/max(U2017res['[cP7]']))
Out[74]:
[<matplotlib.lines.Line2D at 0x40b56773d0>]
In [75]:
%pylab inline
U2017 = te.loada('Models/U2019.4.txt')
#
U2017.reset() #= te.loada('Models/U2017.txt')


U2017.setValue('ge3',1)
range(0,24*10,1)
res = U2017.simulate(0,24*15,24*15)
#plot(res['time'],res['[cLUX]']/max(res['[cLUX]']))
plot(res['time'],res['[cEC]']/max(res['[cEC]']))
U2017 = te.loada('Models/U2019.4.txt')

U2017.setValue('ge3',100)
#xlim(-24,24*5)
xticks(range(0,24*5,24))
U2017.setValue('cE3c',U2017.getValue('cE3c')*U2017.getValue('ge3'))
U2017.setValue('cE3n',U2017.getValue('cE3n')*U2017.getValue('ge3'))
U2017.setValue('cE4',U2017.getValue('cE4')*U2017.getValue('ge3'))
U2017.setValue('cEG',U2017.getValue('cEG')*U2017.getValue('ge3'))
U2017.setValue('cEGn',U2017.getValue('cEGn')*U2017.getValue('ge3'))
U2017.setValue('cLUX',U2017.getValue('cLUX')*U2017.getValue('ge3'))
U2017.setValue('cGc',U2017.getValue('cGc')*U2017.getValue('ge3'))
U2017.setValue('cGn',U2017.getValue('cGn')*U2017.getValue('ge3'))
U2017.setValue('cEC',U2017.getValue('cEC')*U2017.getValue('ge3'))
U2017.setValue('cZTL',U2017.getValue('cZTL')*U2017.getValue('ge3'))
U2017.setValue('cZG',U2017.getValue('cZG')*U2017.getValue('ge3'))
U2017.setValue('cE34',U2017.getValue('cE34')*U2017.getValue('ge3'))
res = U2017.simulate(0,24*15,300)
#plot(res['time'],res['[cLUX]']/max(res['[cLUX]']))
plot(res['time'],res['[cEC]']/max(res['[cEC]']))
#plot(res['time'],res['[cLUX]'])

Populating the interactive namespace from numpy and matplotlib
Out[75]:
[<matplotlib.lines.Line2D at 0x40bc531dd0>]
In [76]:
m = pro_data['cL'].values()[0][0] 
for i in pro_data['cL'].values():
    if i[0] > m:
        m = i[0]
m
Out[76]:
106759.4287181982
In [77]:
U2017.resetAll()
def residuals(params, U2017, id,gene, pro_data):
    U2017.resetAll()
    U2017.setValue(id, params[id].value)
    U2017.setValue(gene,U2017.getValue(gene)*params[id].value)
    res = U2017.simulate(0,24*15,24*15)
    
    max(res['['+gene+']'])
    
    residuals = []
    pred = {}
    for idx,i in enumerate(res['['+gene+']']):
        pred[around(res['time'][idx])-24*5] = res['['+gene+']'][idx]
    for i in pro_data[gene].keys():
            residuals.append(((pro_data[gene][i][0])-(pred[int(i)]))**2)
    return residuals

def residuals_max(params, U2017, id,gene, pro_data):
    U2017.resetAll()
    U2017.setValue(id, params[id].value)
    U2017.setValue(gene,U2017.getValue(gene)*params[id].value)
    res = U2017.simulate(0,24*15,24*15)
    
    max(res['['+gene+']'])
    
    residuals = []
    pred = {}
    for idx,i in enumerate(res['['+gene+']']):
        pred[around(res['time'][idx])-24*5] = res['['+gene+']'][idx]
    
    m1 = pro_data[gene].values()[0][0] 
    for i in pro_data[gene].values():
        if i[0] > m1:
            m1 = i[0]
    
    for i in pro_data[gene].keys():
            residuals.append(((pro_data[gene][i][0])-max(res['['+gene+']']))**2)
    return (m1-max(res['['+gene+']']))**2
In [78]:
p = {}
params = Parameters()
params.add('gl',value=1e5, min=1e4, max=1e6, vary=True)
out = minimize(residuals_max,params,args=(U2017,'gl','cL',pro_data), method='leastsq')
report_fit(out.params)
p['gl'] = out.params['gl'].value

[[Variables]]
    gl:  109962.841 +/- 0.00172892 (0.00%) (init = 100000)
In [79]:
figure(dpi=600)
lh = []
r = linspace(1,0.2e6,100)
for i in r:
    params['gl'].value=i
    lh.append(sum(residuals_max(params, U2017, 'gl','cL', pro_data)))

yticks(fontsize=16)
xticks(fontsize=16)
plot(r,lh)
ylabel('Cost',fontsize=22)
xlabel('scale factor value',fontsize=22)
title('Cost landscape', fontsize=22)
tight_layout()
xscale('log')
savefig('cost_landscape_cL_max_U2019_4.svg', format='svg',dpi=600, transparent=True)
In [80]:
params = Parameters()
params.add('gp7',value=6e4, min=1e4, max=1e6, vary=True)
out = minimize(residuals,params,args=(U2017,'gp7','cP7',pro_data), method='leastsq')
report_fit(out.params)
p['gp7'] = out.params['gp7'].value

[[Variables]]
    gp7:  72532.2121 +/- 140.993082 (0.19%) (init = 60000)
In [81]:
figure(dpi=600)
lh = []
r = linspace(1,0.8e5,100)
for i in r:
    params['gp7'].value=i
    lh.append(sum(residuals(params, U2017, 'gp7','cP7', pro_data)))

plot(r,lh)

yticks(fontsize=16)
xticks(fontsize=16)
plot(r,lh)
ylabel('Cost',fontsize=22)
xlabel('scale factor value',fontsize=22)
title('Cost landscape', fontsize=22)
tight_layout()
xscale('log')
savefig('cost_landscape_P7_max_U2019_4.svg', format='svg',dpi=600, transparent=True)
In [87]:
## this was commented as the residual function was changed for the other genes up, here several genes had to be increased perhaps this has to be adjusted up, the issue is that the general paramters vlues where removed
def residuals(params, U2017, id,gene, pro_data):
    U2017.resetAll()
    #for i in U2017p.keys():
    #    U2017.setValue(i,U2017p.getByKey(i))
    U2017.setValue(id, params[id].value)
    U2017.setValue(gene,U2017.getValue(gene)*params[id].value)
    U2017.setValue('cE3c',U2017.getValue('cE3c')*params[id].value)
    U2017.setValue('cE3n',U2017.getValue('cE3n')*params[id].value)
    U2017.setValue('cE4',U2017.getValue('cE4')*params[id].value)
    U2017.setValue('cEGc',U2017.getValue('cEGc')*params[id].value)
    U2017.setValue('cEGn',U2017.getValue('cEGn')*params[id].value)
    #U2017.setValue('cLUX',U2017.getValue('cLUX')*params[id].value)
    U2017.setValue('cGc',U2017.getValue('cGc')*params[id].value)
    U2017.setValue('cGn',U2017.getValue('cGn')*params[id].value)
    U2017.setValue('cEC',U2017.getValue('cEC')*params[id].value)
    U2017.setValue('cZTL',U2017.getValue('cZTL')*params[id].value)
    U2017.setValue('cZG',U2017.getValue('cZG')*params[id].value)
    U2017.setValue('cE34',U2017.getValue('cE34')*params[id].value)
    res = U2017.simulate(0,24*15,24*15)
    residuals = []
    pred = {}
    for idx,i in enumerate(res['['+gene+']']):
        pred[around(res['time'][idx])-24*5] = res['['+gene+']'][idx]
    for i in pro_data[gene].keys():
            residuals.append(((pro_data[gene][i][0])-(pred[int(i)]))**2)
    return residuals
In [92]:
params = Parameters()
params.add('ge3',value=3e4, min=1e4, max=1e6, vary=True)
out = minimize(residuals,params,args=(U2017,'ge3','cLUX',pro_data), method='leastsq')
report_fit(out.params)
params = out.params
p['ge3'] = out.params['ge3'].value

[[Variables]]
    ge3:  29717.8208 +/- 54.8892322 (0.18%) (init = 30000)
In [90]:
figure(dpi=600)
lh = []
r = linspace(1,1.2e5,100)
for i in r:
    params['ge3'].value=i
    lh.append(sum(residuals(params, U2017, 'ge3','cLUX', pro_data)))

plot(r,lh)

yticks(fontsize=16)
xticks(fontsize=16)
plot(r,lh)
ylabel('Cost',fontsize=22)
xlabel('scale factor value',fontsize=22)
title('Cost landscape', fontsize=22)
tight_layout()
xscale('log')
#savefig('cost_landscape_E3_max_U2019_4_U2019_4.svg', format='svg',dpi=600, transparent=True)
In [101]:
U2017_scaling_factors = DataFrame.from_dict(p, orient="index")
U2017_scaling_factors.to_excel("U2019_4_scaling_factors.xls")
In [102]:
p_loaded = read_excel("U2019_4_scaling_factors.xls")
In [104]:
sf={}
for idx, scaling_factor in enumerate(p_loaded.iloc[:,0]):
    sf[scaling_factor]= p_loaded.iloc[idx,1]
sf
Out[104]:
{u'ge3': 29717.82075330698,
 u'gl': 109962.84089035697,
 u'gp7': 72532.2121212713}
In [105]:
U2017 = te.loada('Models/U2019.4.txt')

for i in p.keys():
    U2017.setValue(i,sf[i])
#U2017.setValue('m15',0.3)
res = U2017.simulate(0,24*15,3000)

Here I create the files for the scaled versions of the model, which in tellurium version and sbml. The sbml version will need to have the special light function

In [130]:
#U2017.exportToAntimony("Models/U2019.4.scaled.txt")
#U2017.exportToSBML("Models/U2019.4.scaled.sbml")
In [117]:
import pickle
In [118]:
with open('pro_data_predicted_from_mRNA.pkl', 'wb') as f:
    pickle.dump(pro_data, f)
        
with open('pro_data_predicted_from_mRNA.pkl', 'rb') as f:
    pro_data_loaded = pickle.load(f)
In [119]:
pro_data_loaded
Out[119]:
{'cE3': {96: (37922.930332050455, 0.1),
  97: (36035.428651557595, 0.1),
  98: (34033.27873077198, 0.1),
  99: (31960.43029922922, 0.1),
  100: (29922.804717238476, 0.1),
  101: (28105.170663355228, 0.1),
  102: (26782.359687776177, 0.1),
  103: (26280.035933815914, 0.1),
  104: (26770.490540796138, 0.1),
  105: (28130.193633629246, 0.1),
  106: (30057.969905236823, 0.1),
  107: (32418.025690299102, 0.1),
  108: (35039.259716459936, 0.1),
  109: (37399.89165126717, 0.1),
  110: (39430.59051256655, 0.1),
  111: (42073.30070006311, 0.1),
  112: (44815.31634904898, 0.1),
  113: (44577.59815691232, 0.1),
  114: (42335.18160338078, 0.1),
  115: (40202.27013180229, 0.1),
  116: (39735.96275370937, 0.1),
  117: (40600.6816895264, 0.1),
  118: (40423.91870780635, 0.1),
  119: (39613.42363996578, 0.1),
  120: (38763.6416395758, 0.1),
  121: (37339.33441543231, 0.1),
  122: (35290.9321113852, 0.1),
  123: (33138.73085174258, 0.1),
  124: (31227.8523819894, 0.1),
  125: (29934.48590166477, 0.1),
  126: (29879.348100537525, 0.1),
  127: (31160.078511537347, 0.1),
  128: (33188.063192087844, 0.1),
  129: (36402.47345707288, 0.1),
  130: (40603.87353032368, 0.1),
  131: (42460.4246623872, 0.1),
  132: (41823.69303135002, 0.1),
  133: (41420.843832462146, 0.1),
  134: (43716.21045587089, 0.1),
  135: (46691.3604304598, 0.1),
  136: (47108.42220792403, 0.1),
  137: (46824.256392899944, 0.1),
  138: (47883.160003135075, 0.1),
  139: (50471.29040510442, 0.1),
  140: (53055.4255143393, 0.1),
  141: (54609.66360213741, 0.1),
  142: (54896.86045114923, 0.1),
  143: (54059.27708344895, 0.1),
  144: (52671.13810827745, 0.1),
  145: (51649.89791811762, 0.1),
  146: (51193.12671857809, 0.1),
  147: (49749.283001770746, 0.1),
  148: (47248.05679089695, 0.1),
  149: (45012.33912301276, 0.1),
  150: (43913.552243367136, 0.1),
  151: (42413.747358379995, 0.1),
  152: (39982.104083784114, 0.1),
  153: (37962.75637713644, 0.1),
  154: (38983.66781536405, 0.1),
  155: (44084.054602632255, 0.1),
  156: (47582.974085704926, 0.1),
  157: (49295.56035285918, 0.1),
  158: (50985.262710952, 0.1),
  159: (52554.516034824024, 0.1),
  160: (53430.136345445, 0.1),
  161: (53839.989482521174, 0.1),
  162: (54267.85012590509, 0.1),
  163: (55048.53890355082, 0.1),
  164: (56220.810775900456, 0.1),
  165: (57307.21115377372, 0.1),
  166: (57426.7261420727, 0.1),
  167: (55914.43993753194, 0.1)},
 'cE4': {96: (4566.5545671258005, 0.1),
  97: (3271.0572539643094, 0.1),
  98: (2544.622424024666, 0.1),
  99: (2111.2613709740767, 0.1),
  100: (1697.2841705841463, 0.1),
  101: (1436.555840041637, 0.1),
  102: (1726.302352167201, 0.1),
  103: (4653.290792398484, 0.1),
  104: (18361.607314504297, 0.1),
  105: (45915.36312788802, 0.1),
  106: (70980.32012860113, 0.1),
  107: (92285.99669196617, 0.1),
  108: (107514.47154546333, 0.1),
  109: (98681.29711370713, 0.1),
  110: (76702.08837669816, 0.1),
  111: (57782.6083767403, 0.1),
  112: (44086.51902583135, 0.1),
  113: (32253.193551807006, 0.1),
  114: (22617.535198773854, 0.1),
  115: (15908.260073310688, 0.1),
  116: (11979.548655746985, 0.1),
  117: (10169.298770854399, 0.1),
  118: (8538.356699259153, 0.1),
  119: (6422.549275609669, 0.1),
  120: (4538.88996526022, 0.1),
  121: (3171.964802908097, 0.1),
  122: (2238.385318541816, 0.1),
  123: (1589.9791627351376, 0.1),
  124: (1130.7992720631244, 0.1),
  125: (838.7992892831464, 0.1),
  126: (815.9499891144026, 0.1),
  127: (2025.9626381794578, 0.1),
  128: (10392.362728105654, 0.1),
  129: (41300.68155496551, 0.1),
  130: (91196.41082837667, 0.1),
  131: (115602.67109012198, 0.1),
  132: (107096.2862059866, 0.1),
  133: (91329.30006358733, 0.1),
  134: (79742.05550541102, 0.1),
  135: (71155.65308377414, 0.1),
  136: (61804.12151821624, 0.1),
  137: (50862.97566834542, 0.1),
  138: (40257.7102814985, 0.1),
  139: (32048.65286235233, 0.1),
  140: (26364.237664813467, 0.1),
  141: (21540.694711491073, 0.1),
  142: (16866.92598139636, 0.1),
  143: (12894.751341629724, 0.1),
  144: (9855.217405184183, 0.1),
  145: (7598.605352090931, 0.1),
  146: (5870.679674441468, 0.1),
  147: (4458.40737828358, 0.1),
  148: (3350.630865233303, 0.1),
  149: (2645.758012331802, 0.1),
  150: (2396.9663005054745, 0.1),
  151: (2332.6027294481632, 0.1),
  152: (2346.0298022860607, 0.1),
  153: (3722.983688456442, 0.1),
  154: (12688.257850824299, 0.1),
  155: (30886.730186785975, 0.1),
  156: (42099.344743685746, 0.1),
  157: (47326.10340629515, 0.1),
  158: (51397.83391818968, 0.1),
  159: (50180.43782731251, 0.1),
  160: (43871.047766385236, 0.1),
  161: (38184.60654764137, 0.1),
  162: (35767.84092704913, 0.1),
  163: (34301.85678818751, 0.1),
  164: (30901.125790550213, 0.1),
  165: (25953.2376959809, 0.1),
  166: (20938.112092356587, 0.1),
  167: (16916.67802166559, 0.1)},
 'cL': {96: (58910.554068438825, 0.1),
  97: (89506.78137150322, 0.1),
  98: (86570.0115883527, 0.1),
  99: (29303.653739416193, 0.1),
  100: (6756.7643117377, 0.1),
  101: (4508.483065710765, 0.1),
  102: (17060.74433520031, 0.1),
  103: (11763.118667172019, 0.1),
  104: (2762.9003697488315, 0.1),
  105: (734.3273254214948, 0.1),
  106: (629.9735975527343, 0.1),
  107: (688.210657694042, 0.1),
  108: (490.81493034538875, 0.1),
  109: (473.36300956419655, 0.1),
  110: (693.2136894467521, 0.1),
  111: (735.4596503114919, 0.1),
  112: (542.9313786121231, 0.1),
  113: (458.08980568924085, 0.1),
  114: (756.7741783114504, 0.1),
  115: (3255.2981951130455, 0.1),
  116: (16410.552203697902, 0.1),
  117: (34914.752401614285, 0.1),
  118: (42404.10531174643, 0.1),
  119: (67591.93261641228, 0.1),
  120: (106759.4287181982, 0.1),
  121: (70339.46919568657, 0.1),
  122: (28343.49167939349, 0.1),
  123: (11295.41563746694, 0.1),
  124: (5457.272513319294, 0.1),
  125: (4980.51496639676, 0.1),
  126: (6939.2858139126965, 0.1),
  127: (4848.913455604003, 0.1),
  128: (1860.9015571999323, 0.1),
  129: (826.4485149143803, 0.1),
  130: (509.9213285409796, 0.1),
  131: (212.54339785082692, 0.1),
  132: (71.66115507231221, 0.1),
  133: (68.63294999692496, 0.1),
  134: (264.20637310937207, 0.1),
  135: (480.5986576990896, 0.1),
  136: (449.27168415340964, 0.1),
  137: (757.4342644315112, 0.1),
  138: (2941.991676986211, 0.1),
  139: (9204.747462128018, 0.1),
  140: (16915.255830657454, 0.1),
  141: (21508.212642588544, 0.1),
  142: (24638.088093966286, 0.1),
  143: (31572.95413358934, 0.1),
  144: (46416.35289658086, 0.1),
  145: (58786.64471703827, 0.1),
  146: (58158.854609986796, 0.1),
  147: (56143.6503018239, 0.1),
  148: (49851.64943579772, 0.1),
  149: (33384.823066093755, 0.1),
  150: (21457.514217600787, 0.1),
  151: (24950.75588268126, 0.1),
  152: (35746.36145699529, 0.1),
  153: (21400.869812201694, 0.1),
  154: (6440.545095128797, 0.1),
  155: (1728.9135479278955, 0.1),
  156: (626.0205012851876, 0.1),
  157: (449.6812961669427, 0.1),
  158: (603.2623788017484, 0.1),
  159: (909.6309886689648, 0.1),
  160: (1371.7540315746296, 0.1),
  161: (2501.4994905467747, 0.1),
  162: (5394.978703272053, 0.1),
  163: (11119.185884439263, 0.1),
  164: (19545.576200767682, 0.1),
  165: (28421.947478703773, 0.1),
  166: (33336.94470480675, 0.1),
  167: (30756.158493848416, 0.1)},
 'cLUX': {96: (70649.94604933374, 0.1),
  97: (68403.92170012437, 0.1),
  98: (66178.52141009286, 0.1),
  99: (63977.54112552265, 0.1),
  100: (61829.86885217761, 0.1),
  101: (59975.762988224626, 0.1),
  102: (59001.80089699478, 0.1),
  103: (59936.15075344081, 0.1),
  104: (63953.60521997696, 0.1),
  105: (71752.10111323271, 0.1),
  106: (82309.24772230642, 0.1),
  107: (92289.83904186949, 0.1),
  108: (98192.875631402, 0.1),
  109: (99230.63653348955, 0.1),
  110: (97323.56537844845, 0.1),
  111: (94321.21029062626, 0.1),
  112: (90906.16854189664, 0.1),
  113: (87314.08625754385, 0.1),
  114: (83769.0077067292, 0.1),
  115: (80380.76612245668, 0.1),
  116: (77337.15153227691, 0.1),
  117: (75153.01740683417, 0.1),
  118: (73809.42111158071, 0.1),
  119: (72273.03533385218, 0.1),
  120: (70135.0818106056, 0.1),
  121: (67742.24195774952, 0.1),
  122: (65356.54687233729, 0.1),
  123: (63097.758116011224, 0.1),
  124: (61034.26223670364, 0.1),
  125: (59179.614137073564, 0.1),
  126: (57675.647611735796, 0.1),
  127: (57535.76597543213, 0.1),
  128: (62543.185265790664, 0.1),
  129: (77558.65212043938, 0.1),
  130: (96939.9992726125, 0.1),
  131: (108419.59499932155, 0.1),
  132: (111433.32975802754, 0.1),
  133: (111452.84587023094, 0.1),
  134: (111497.31868485887, 0.1),
  135: (112969.39739504445, 0.1),
  136: (115689.99863065746, 0.1),
  137: (116841.44072113249, 0.1),
  138: (115511.65084305947, 0.1),
  139: (113409.02262869288, 0.1),
  140: (111722.81337717257, 0.1),
  141: (110464.07473148129, 0.1),
  142: (108865.31122066107, 0.1),
  143: (106364.46947617644, 0.1),
  144: (103272.88519687524, 0.1),
  145: (100351.24309621811, 0.1),
  146: (98101.10694738278, 0.1),
  147: (95881.41515550723, 0.1),
  148: (93195.03433158652, 0.1),
  149: (90735.27061503354, 0.1),
  150: (89379.74504351914, 0.1),
  151: (88749.81016709507, 0.1),
  152: (88048.14539936582, 0.1),
  153: (88928.40985501274, 0.1),
  154: (95363.28395824981, 0.1),
  155: (106528.23756302256, 0.1),
  156: (114980.38845824832, 0.1),
  157: (119369.24048392124, 0.1),
  158: (121246.30176595795, 0.1),
  159: (121469.80178563463, 0.1),
  160: (120829.81396235638, 0.1),
  161: (120406.33129398199, 0.1),
  162: (120974.40300529079, 0.1),
  163: (122020.29904493649, 0.1),
  164: (122351.32696722876, 0.1),
  165: (121455.01765551591, 0.1),
  166: (119481.19470811926, 0.1),
  167: (116835.91272031, 0.1)},
 'cP5': {96: (2687.5327529136475, 0.1),
  97: (2223.488190257279, 0.1),
  98: (2058.6517468454167, 0.1),
  99: (2088.9578016022147, 0.1),
  100: (2276.91074602263, 0.1),
  101: (2691.7538147704518, 0.1),
  102: (3809.6459626412534, 0.1),
  103: (8294.123488492123, 0.1),
  104: (23474.88081424806, 0.1),
  105: (41948.42991531756, 0.1),
  106: (47206.54798009061, 0.1),
  107: (44056.92391633113, 0.1),
  108: (36419.394488452366, 0.1),
  109: (26217.62139599148, 0.1),
  110: (18086.21200933459, 0.1),
  111: (12312.57187715425, 0.1),
  112: (8122.755677172437, 0.1),
  113: (5330.022080896432, 0.1),
  114: (3605.2028403936934, 0.1),
  115: (2584.835087680827, 0.1),
  116: (2085.9965921093153, 0.1),
  117: (2338.568309545493, 0.1),
  118: (3476.4697239460584, 0.1),
  119: (3605.766023285633, 0.1),
  120: (2860.987181550808, 0.1),
  121: (2642.087435731194, 0.1),
  122: (2436.5318608701377, 0.1),
  123: (2243.6662589266025, 0.1),
  124: (2053.631154446891, 0.1),
  125: (1984.8775389728505, 0.1),
  126: (2577.3089652144745, 0.1),
  127: (7177.338559136989, 0.1),
  128: (26597.399972525956, 0.1),
  129: (52372.25021141758, 0.1),
  130: (60816.64794315633, 0.1),
  131: (58768.818903026186, 0.1),
  132: (53807.44026667461, 0.1),
  133: (48193.49873389751, 0.1),
  134: (42903.47334270371, 0.1),
  135: (38808.36236610994, 0.1),
  136: (36353.93924757793, 0.1),
  137: (34239.49959964348, 0.1),
  138: (31373.94975488764, 0.1),
  139: (28411.97448374716, 0.1),
  140: (25918.764229110056, 0.1),
  141: (23804.7678421133, 0.1),
  142: (21595.005746949417, 0.1),
  143: (19024.806255718566, 0.1),
  144: (16443.811571730377, 0.1),
  145: (14218.280077966272, 0.1),
  146: (12494.060611686942, 0.1),
  147: (11080.444408866386, 0.1),
  148: (9811.068861353266, 0.1),
  149: (9052.865536943014, 0.1),
  150: (9527.569229749493, 0.1),
  151: (10019.360613805204, 0.1),
  152: (9671.093205716601, 0.1),
  153: (11005.883914190208, 0.1),
  154: (23011.83120247067, 0.1),
  155: (34468.77076931919, 0.1),
  156: (32927.856365269196, 0.1),
  157: (29477.45040619848, 0.1),
  158: (27455.14107356111, 0.1),
  159: (26666.442291402203, 0.1),
  160: (25419.629957602094, 0.1),
  161: (23924.97475173252, 0.1),
  162: (22777.40132574089, 0.1),
  163: (21893.297304549822, 0.1),
  164: (20899.40275294496, 0.1),
  165: (19599.188454199142, 0.1),
  166: (17997.377199693, 0.1),
  167: (16210.123285547665, 0.1)},
 'cP7': {96: (1719.9822193713842, 0.1),
  97: (1613.6215581210215, 0.1),
  98: (1701.4641625165875, 0.1),
  99: (2362.805823643368, 0.1),
  100: (4773.039648681314, 0.1),
  101: (8160.055388453035, 0.1),
  102: (10432.962564588417, 0.1),
  103: (13597.260732694836, 0.1),
  104: (21772.785936822507, 0.1),
  105: (32531.44906502603, 0.1),
  106: (37342.956561959705, 0.1),
  107: (36223.690853744214, 0.1),
  108: (31936.032751107843, 0.1),
  109: (26200.15356783291, 0.1),
  110: (20748.025692203464, 0.1),
  111: (16243.104037135792, 0.1),
  112: (12624.452037071527, 0.1),
  113: (9724.512327435586, 0.1),
  114: (7449.125009783082, 0.1),
  115: (5701.044271426131, 0.1),
  116: (4372.75333090303, 0.1),
  117: (3366.749243027449, 0.1),
  118: (2604.9453481397036, 0.1),
  119: (2044.026264898816, 0.1),
  120: (1707.473534117763, 0.1),
  121: (1971.170152813368, 0.1),
  122: (5017.880827468266, 0.1),
  123: (14391.560037536292, 0.1),
  124: (24770.69011150918, 0.1),
  125: (29887.330691658713, 0.1),
  126: (31423.740005984277, 0.1),
  127: (34680.011645612634, 0.1),
  128: (45214.02939906451, 0.1),
  129: (56834.62708248042, 0.1),
  130: (59200.41446130011, 0.1),
  131: (55518.94391431924, 0.1),
  132: (49622.33982444178, 0.1),
  133: (42953.48446961325, 0.1),
  134: (36620.277221607634, 0.1),
  135: (31390.71696488101, 0.1),
  136: (27526.195388675707, 0.1),
  137: (24045.851686677193, 0.1),
  138: (20479.441416254456, 0.1),
  139: (17287.086548441057, 0.1),
  140: (14607.042866312948, 0.1),
  141: (12375.939747608278, 0.1),
  142: (10525.575094870059, 0.1),
  143: (9093.6278941361, 0.1),
  144: (8199.89637329123, 0.1),
  145: (7758.1608967336015, 0.1),
  146: (7605.10186310404, 0.1),
  147: (8367.455179491977, 0.1),
  148: (12082.183394609303, 0.1),
  149: (20788.36939329959, 0.1),
  150: (32145.23210692608, 0.1),
  151: (40607.240927497245, 0.1),
  152: (44979.46264721154, 0.1),
  153: (49454.09143137253, 0.1),
  154: (56839.0839819828, 0.1),
  155: (61713.08550720408, 0.1),
  156: (59573.345930997704, 0.1),
  157: (53425.287784487875, 0.1),
  158: (46359.22528499026, 0.1),
  159: (39716.34095262406, 0.1),
  160: (33895.96234933267, 0.1),
  161: (28922.294188974607, 0.1),
  162: (24689.6181072745, 0.1),
  163: (21047.91471826665, 0.1),
  164: (17905.657516770774, 0.1),
  165: (15247.334952729652, 0.1),
  166: (13117.316619830812, 0.1),
  167: (11836.525165491536, 0.1)},
 'cP9': {96: (288.3247258782284, 0.1),
  97: (214.5560648073497, 0.1),
  98: (1521.2802369804363, 0.1),
  99: (3124.1718618044365, 0.1),
  100: (2391.67173460836, 0.1),
  101: (3599.020795599749, 0.1),
  102: (7240.359034919757, 0.1),
  103: (3779.8127533836314, 0.1),
  104: (858.4834247598308, 0.1),
  105: (285.66024809206476, 0.1),
  106: (299.98700949149463, 0.1),
  107: (288.3308484734385, 0.1),
  108: (204.77714495180788, 0.1),
  109: (221.485066630427, 0.1),
  110: (235.20228365072225, 0.1),
  111: (89.44594413923201, 0.1),
  112: (21.57193804172612, 0.1),
  113: (18.3488368720072, 0.1),
  114: (50.134804850321636, 0.1),
  115: (29.618616680977436, 0.1),
  116: (10.951683132688295, 0.1),
  117: (108.00158762739343, 0.1),
  118: (1660.0019017512518, 0.1),
  119: (562.0544983766016, 0.1),
  120: (64.94243770939711, 0.1),
  121: (17.760765438288317, 0.1),
  122: (1342.7591245291628, 0.1),
  123: (8309.933534915614, 0.1),
  124: (5442.230565661302, 0.1),
  125: (3718.98490741422, 0.1),
  126: (4411.588826978921, 0.1),
  127: (2716.188466423511, 0.1),
  128: (800.4828339674668, 0.1),
  129: (188.0015047076733, 0.1),
  130: (49.16322017398147, 0.1),
  131: (17.386071578984865, 0.1),
  132: (8.742293928590394, 0.1),
  133: (5.704823571948542, 0.1),
  134: (5.346221152727214, 0.1),
  135: (9.022000821369831, 0.1),
  136: (13.94123543232016, 0.1),
  137: (8.040209500486737, 0.1),
  138: (2.7953553041414025, 0.1),
  139: (2.1435299315701104, 0.1),
  140: (4.952252006198682, 0.1),
  141: (8.47607312566374, 0.1),
  142: (9.20863503656202, 0.1),
  143: (11.308460642685395, 0.1),
  144: (25.05188001527406, 0.1),
  145: (119.1152627080142, 0.1),
  146: (572.8994159321272, 0.1),
  147: (1091.673807036661, 0.1),
  148: (1042.936476174327, 0.1),
  149: (1056.3346693735969, 0.1),
  150: (1643.8385012095428, 0.1),
  151: (3566.466137619713, 0.1),
  152: (5759.562854637112, 0.1),
  153: (3617.5986704255265, 0.1),
  154: (1082.1265506676755, 0.1),
  155: (256.143086471718, 0.1),
  156: (73.99800759589107, 0.1),
  157: (43.236380963341894, 0.1),
  158: (40.14634880498924, 0.1),
  159: (23.854760645868264, 0.1),
  160: (9.813945742711146, 0.1),
  161: (5.648486294626583, 0.1),
  162: (6.027882265219004, 0.1),
  163: (7.431029383212584, 0.1),
  164: (8.403563092316821, 0.1),
  165: (9.190268439684132, 0.1),
  166: (10.634523504576775, 0.1),
  167: (14.348849232105506, 0.1)},
 'cT': {96: (2093.105500592261, 0.1),
  97: (1396.2573958724738, 0.1),
  98: (1008.0583405251457, 0.1),
  99: (701.8881990508224, 0.1),
  100: (493.6775524030296, 0.1),
  101: (499.3327067472529, 0.1),
  102: (958.2375236314255, 0.1),
  103: (2642.2950875375664, 0.1),
  104: (7066.09137187075, 0.1),
  105: (14235.256706198632, 0.1),
  106: (20686.553304645928, 0.1),
  107: (23313.323802679766, 0.1),
  108: (21489.00000175553, 0.1),
  109: (16719.778998179212, 0.1),
  110: (11762.203371399948, 0.1),
  111: (8292.885920017805, 0.1),
  112: (5894.4628990561005, 0.1),
  113: (3791.5706645040364, 0.1),
  114: (2403.743322581635, 0.1),
  115: (2264.0704107522365, 0.1),
  116: (3651.1007077551276, 0.1),
  117: (5045.618971278726, 0.1),
  118: (4612.396778077561, 0.1),
  119: (3310.139671736735, 0.1),
  120: (2069.4496094983083, 0.1),
  121: (1144.8639711015642, 0.1),
  122: (622.2065610199459, 0.1),
  123: (452.94743129192386, 0.1),
  124: (716.9318379649814, 0.1),
  125: (2497.486310323754, 0.1),
  126: (5679.728505171122, 0.1),
  127: (3625.80010709758, 0.1),
  128: (1427.48397574397, 0.1),
  129: (1068.6925766502427, 0.1),
  130: (9549.453540861357, 0.1),
  131: (29981.58067324995, 0.1),
  132: (25406.956760182235, 0.1),
  133: (16625.788213865642, 0.1),
  134: (13445.794428079085, 0.1),
  135: (12505.368368272088, 0.1),
  136: (11011.175468020852, 0.1),
  137: (8968.041922238008, 0.1),
  138: (7297.9632631186, 0.1),
  139: (6744.490306919582, 0.1),
  140: (6885.468783815842, 0.1),
  141: (6229.075922507217, 0.1),
  142: (4895.482854527197, 0.1),
  143: (4063.826344572895, 0.1),
  144: (3863.6351886855496, 0.1),
  145: (3601.763523096394, 0.1),
  146: (2985.946015130843, 0.1),
  147: (2301.6707695763553, 0.1),
  148: (1877.9611559336722, 0.1),
  149: (1983.4884249942124, 0.1),
  150: (2602.0839003051433, 0.1),
  151: (2991.0317964240267, 0.1),
  152: (3120.3559756642967, 0.1),
  153: (4555.484185711643, 0.1),
  154: (9682.333958240904, 0.1),
  155: (17155.81821513042, 0.1),
  156: (21093.347384181823, 0.1),
  157: (20565.498013971985, 0.1),
  158: (17360.75016158397, 0.1),
  159: (13020.762534320895, 0.1),
  160: (9520.55130137989, 0.1),
  161: (8271.094701588792, 0.1),
  162: (9912.040001223886, 0.1),
  163: (15344.344292767957, 0.1),
  164: (23062.421775594525, 0.1),
  165: (24021.093754308724, 0.1),
  166: (14372.63873888575, 0.1),
  167: (5462.68067156176, 0.1)}}

Conversion of the protein data into a excel file for further use by others. The pickle file can also be read using python however it woule be difficult to access by futher users

In [202]:
#This transforms the predicted data into DataFrame that can then be sevaed as an excel file and uploaded
#into FAIRDOM hub
pro_data_matrix = zeros((len(pro_data["cL"]),len(pro_data.keys())+1))
header=[]
header.append("time")
for idx, key in enumerate(pro_data.keys()):
    header.append(key)
    for tp in sort(pro_data[key].keys()):
        pro_data_matrix[tp-96,0]=tp
        pro_data_matrix[tp-96,idx+1]=pro_data[key][tp][0]
pro_data_dataframe_predictions = DataFrame(data=pro_data_matrix,columns=header) 
pro_data_dataframe_predictions.to_excel("Urquiza2016_protein_predictions_from_TiMet.xls")
In [197]:
for i in range(1,9):
    plot(pro_data_matrix[:,0],pro_data_matrix[:,i], label=header[i])
legend(loc="upper right")
Out[197]:
<matplotlib.legend.Legend at 0x40fc116fd0>
In [199]:
Out[199]:
time cL cP9 cLUX cE4 cE3 cT cP5 cP7
0 96.0 58910.554068 288.324726 70649.946049 4566.554567 37922.930332 2093.105501 2687.532753 1719.982219
1 97.0 89506.781372 214.556065 68403.921700 3271.057254 36035.428652 1396.257396 2223.488190 1613.621558
2 98.0 86570.011588 1521.280237 66178.521410 2544.622424 34033.278731 1008.058341 2058.651747 1701.464163
3 99.0 29303.653739 3124.171862 63977.541126 2111.261371 31960.430299 701.888199 2088.957802 2362.805824
4 100.0 6756.764312 2391.671735 61829.868852 1697.284171 29922.804717 493.677552 2276.910746 4773.039649
5 101.0 4508.483066 3599.020796 59975.762988 1436.555840 28105.170663 499.332707 2691.753815 8160.055388
6 102.0 17060.744335 7240.359035 59001.800897 1726.302352 26782.359688 958.237524 3809.645963 10432.962565
7 103.0 11763.118667 3779.812753 59936.150753 4653.290792 26280.035934 2642.295088 8294.123488 13597.260733
8 104.0 2762.900370 858.483425 63953.605220 18361.607315 26770.490541 7066.091372 23474.880814 21772.785937
9 105.0 734.327325 285.660248 71752.101113 45915.363128 28130.193634 14235.256706 41948.429915 32531.449065
10 106.0 629.973598 299.987009 82309.247722 70980.320129 30057.969905 20686.553305 47206.547980 37342.956562
11 107.0 688.210658 288.330848 92289.839042 92285.996692 32418.025690 23313.323803 44056.923916 36223.690854
12 108.0 490.814930 204.777145 98192.875631 107514.471545 35039.259716 21489.000002 36419.394488 31936.032751
13 109.0 473.363010 221.485067 99230.636533 98681.297114 37399.891651 16719.778998 26217.621396 26200.153568
14 110.0 693.213689 235.202284 97323.565378 76702.088377 39430.590513 11762.203371 18086.212009 20748.025692
15 111.0 735.459650 89.445944 94321.210291 57782.608377 42073.300700 8292.885920 12312.571877 16243.104037
16 112.0 542.931379 21.571938 90906.168542 44086.519026 44815.316349 5894.462899 8122.755677 12624.452037
17 113.0 458.089806 18.348837 87314.086258 32253.193552 44577.598157 3791.570665 5330.022081 9724.512327
18 114.0 756.774178 50.134805 83769.007707 22617.535199 42335.181603 2403.743323 3605.202840 7449.125010
19 115.0 3255.298195 29.618617 80380.766122 15908.260073 40202.270132 2264.070411 2584.835088 5701.044271
20 116.0 16410.552204 10.951683 77337.151532 11979.548656 39735.962754 3651.100708 2085.996592 4372.753331
21 117.0 34914.752402 108.001588 75153.017407 10169.298771 40600.681690 5045.618971 2338.568310 3366.749243
22 118.0 42404.105312 1660.001902 73809.421112 8538.356699 40423.918708 4612.396778 3476.469724 2604.945348
23 119.0 67591.932616 562.054498 72273.035334 6422.549276 39613.423640 3310.139672 3605.766023 2044.026265
24 120.0 106759.428718 64.942438 70135.081811 4538.889965 38763.641640 2069.449609 2860.987182 1707.473534
25 121.0 70339.469196 17.760765 67742.241958 3171.964803 37339.334415 1144.863971 2642.087436 1971.170153
26 122.0 28343.491679 1342.759125 65356.546872 2238.385319 35290.932111 622.206561 2436.531861 5017.880827
27 123.0 11295.415637 8309.933535 63097.758116 1589.979163 33138.730852 452.947431 2243.666259 14391.560038
28 124.0 5457.272513 5442.230566 61034.262237 1130.799272 31227.852382 716.931838 2053.631154 24770.690112
29 125.0 4980.514966 3718.984907 59179.614137 838.799289 29934.485902 2497.486310 1984.877539 29887.330692
... ... ... ... ... ... ... ... ... ...
42 138.0 2941.991677 2.795355 115511.650843 40257.710281 47883.160003 7297.963263 31373.949755 20479.441416
43 139.0 9204.747462 2.143530 113409.022629 32048.652862 50471.290405 6744.490307 28411.974484 17287.086548
44 140.0 16915.255831 4.952252 111722.813377 26364.237665 53055.425514 6885.468784 25918.764229 14607.042866
45 141.0 21508.212643 8.476073 110464.074731 21540.694711 54609.663602 6229.075923 23804.767842 12375.939748
46 142.0 24638.088094 9.208635 108865.311221 16866.925981 54896.860451 4895.482855 21595.005747 10525.575095
47 143.0 31572.954134 11.308461 106364.469476 12894.751342 54059.277083 4063.826345 19024.806256 9093.627894
48 144.0 46416.352897 25.051880 103272.885197 9855.217405 52671.138108 3863.635189 16443.811572 8199.896373
49 145.0 58786.644717 119.115263 100351.243096 7598.605352 51649.897918 3601.763523 14218.280078 7758.160897
50 146.0 58158.854610 572.899416 98101.106947 5870.679674 51193.126719 2985.946015 12494.060612 7605.101863
51 147.0 56143.650302 1091.673807 95881.415156 4458.407378 49749.283002 2301.670770 11080.444409 8367.455179
52 148.0 49851.649436 1042.936476 93195.034332 3350.630865 47248.056791 1877.961156 9811.068861 12082.183395
53 149.0 33384.823066 1056.334669 90735.270615 2645.758012 45012.339123 1983.488425 9052.865537 20788.369393
54 150.0 21457.514218 1643.838501 89379.745044 2396.966301 43913.552243 2602.083900 9527.569230 32145.232107
55 151.0 24950.755883 3566.466138 88749.810167 2332.602729 42413.747358 2991.031796 10019.360614 40607.240927
56 152.0 35746.361457 5759.562855 88048.145399 2346.029802 39982.104084 3120.355976 9671.093206 44979.462647
57 153.0 21400.869812 3617.598670 88928.409855 3722.983688 37962.756377 4555.484186 11005.883914 49454.091431
58 154.0 6440.545095 1082.126551 95363.283958 12688.257851 38983.667815 9682.333958 23011.831202 56839.083982
59 155.0 1728.913548 256.143086 106528.237563 30886.730187 44084.054603 17155.818215 34468.770769 61713.085507
60 156.0 626.020501 73.998008 114980.388458 42099.344744 47582.974086 21093.347384 32927.856365 59573.345931
61 157.0 449.681296 43.236381 119369.240484 47326.103406 49295.560353 20565.498014 29477.450406 53425.287784
62 158.0 603.262379 40.146349 121246.301766 51397.833918 50985.262711 17360.750162 27455.141074 46359.225285
63 159.0 909.630989 23.854761 121469.801786 50180.437827 52554.516035 13020.762534 26666.442291 39716.340953
64 160.0 1371.754032 9.813946 120829.813962 43871.047766 53430.136345 9520.551301 25419.629958 33895.962349
65 161.0 2501.499491 5.648486 120406.331294 38184.606548 53839.989483 8271.094702 23924.974752 28922.294189
66 162.0 5394.978703 6.027882 120974.403005 35767.840927 54267.850126 9912.040001 22777.401326 24689.618107
67 163.0 11119.185884 7.431029 122020.299045 34301.856788 55048.538904 15344.344293 21893.297305 21047.914718
68 164.0 19545.576201 8.403563 122351.326967 30901.125791 56220.810776 23062.421776 20899.402753 17905.657517
69 165.0 28421.947479 9.190268 121455.017656 25953.237696 57307.211154 24021.093754 19599.188454 15247.334953
70 166.0 33336.944705 10.634524 119481.194708 20938.112092 57426.726142 14372.638739 17997.377200 13117.316620
71 167.0 30756.158494 14.348849 116835.912720 16916.678022 55914.439938 5462.680672 16210.123286 11836.525165

72 rows × 9 columns

In [111]:
get_max = []
figure(dpi=600)
plot(res['time']-24*11,(res['[cL]']), color='k')
xlim(-24,24*2)
yticks(fontsize=16)
xticks(range(0,24*2,12), fontsize=16)
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(24+12,24*2, color='gray',alpha=0.3)
axvspan(24*2+12,24*3, color='gray',alpha=0.3)
axvspan(12,24, color='gray',alpha=0.3)
d = []
t = []
for i in pro_data['cL'].keys():
    get_max.append(pro_data['cL'][i][0])
    d.append(pro_data['cL'][i][0])
    t.append(i-24*5)
    
d = array(d)
t = array(t)
t = concatenate((t[range(41,72)],t[range(0,40)]))
d = concatenate((d[range(41,72)],d[range(0,40)]))

plot(t,d,'-', color='orange', lw=3)

ylabel('#monomers/cell', fontsize=22)
xlabel('Time hours in constant light (h)', fontsize=22)
yscale('log')
ylim(1e1,1e6)
tight_layout()

savefig('U2019_4_cL_rescaled_total.svg',format='svg',dpi=600, transparent=True)
savefig('U2019_4_cL_rescaled_total.pdf',format='pdf',dpi=600, transparent=True)
print max(get_max[65:])
print min(get_max[65:])
ratios={}
ratios["CL"]=(min(get_max[65:]),max(get_max[65:]))

70339.46919568657
4848.913455604003
In [120]:
def plot_dict(data,var,sf=1):
    time = array(sorted(data[var].keys()))
    d = []
    for i in sorted(data[var].keys()):
        d.append(data[var][i])

    plot(time-24*5,array(d)*sf,color='blue')
In [121]:
figure(dpi=600)
plot(res['time']-24*11,(res['[cP9]']), 'k')
d = []
t = []
for i in pro_data['cP9'].keys():
    #plot(i-24*5,pro_data['cP9'][i][0],'o', color='orange')
    d.append(pro_data['cP9'][i][0])
    t.append(i-24*5)

d = array(d)
t = array(t)
t = concatenate((t[range(41,72)],t[range(0,40)]))
d = concatenate((d[range(41,72)],d[range(0,40)]))    

plot(t,d,'-', color='orange', lw=3)   
#plot_dict(w_data,'cP9',sf=10000)
yticks(fontsize=16)
xlim(-24,24*2)
ylim(1,1e6)
xticks(range(0,24*2,12),fontsize=16)
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(24+12,24*2, color='gray',alpha=0.3)
axvspan(12,24, color='gray',alpha=0.3)
xlabel('Time in constant light (h)', fontsize=22)
ylabel('#monomers/cell', fontsize=22)
yscale('log')
tight_layout()

savefig('U2019_4_cP9_rescaled_total.svg',format='svg',dpi=600, transparent=True)
savefig('U2019_4_cP9_rescaled_total.pdf',format='pdf',dpi=600, transparent=True)
In [122]:
figure(dpi=600)
plot(res['time']-24*11,(res['[cP7]']), 'k')
d = []
t = []
for i in pro_data['cP7'].keys():
    #plot(i-24*5,pro_data['cP9'][i][0],'o', color='orange')
    d.append(pro_data['cP7'][i][0])
    t.append(i-24*5)

d = array(d)
t = array(t)
t = concatenate((t[range(41,72)],t[range(0,40)]))
d = concatenate((d[range(41,72)],d[range(0,40)]))    

plot(t,d,'-', color='orange', lw=3)   
#plot_dict(w_data,'cP9',sf=10000)
yticks(fontsize=16)
xlim(-24,24*2)
ylim(1e2,1e5)
xticks(range(0,24*2,12),fontsize=16)
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(24+12,24*2, color='gray',alpha=0.3)
axvspan(12,24, color='gray',alpha=0.3)
xlabel('Time in constant light (h)', fontsize=22)
ylabel('#monomers/cell', fontsize=22)
yscale('log')
tight_layout()
savefig('U2019_4_cP7_rescaled_total.svg',format='svg',dpi=600, transparent=True)
savefig('U2019_4_cP7_rescaled_total.pdf',format='pdf',dpi=600, transparent=True)
In [123]:
figure(dpi=600)
plot(res['time']-24*11,(res['[cP5]']), 'k')
d = []
t = []
for i in pro_data['cP5'].keys():
    #plot(i-24*5,pro_data['cP9'][i][0],'o', color='orange')
    d.append(pro_data['cP5'][i][0])
    t.append(i-24*5)

d = array(d)
t = array(t)
t = concatenate((t[range(41,72)],t[range(0,40)]))
d = concatenate((d[range(41,72)],d[range(0,40)]))    

plot(t,d,'-', color='orange', lw=3)   
#plot_dict(w_data,'cP9',sf=10000)
yticks(fontsize=16)
xlim(-24,24*2)
ylim(1e2,1e5)
xticks(range(0,24*2,12),fontsize=16)
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(24+12,24*2, color='gray',alpha=0.3)
axvspan(12,24, color='gray',alpha=0.3)
xlabel('Time in constant light (h)', fontsize=22)
ylabel('#monomers/cell', fontsize=22)
yscale('log')
tight_layout()
savefig('U2019_4_cP5_rescaled_total.svg',format='svg',dpi=600, transparent=True)
savefig('U2019_4_cP5_rescaled_total.pdf',format='pdf',dpi=600, transparent=True)
In [124]:
figure(dpi=600)
plot(res['time']-24*11,(res['[cT]']), 'k')
d = []
t = []
for i in pro_data['cT'].keys():
    #plot(i-24*5,pro_data['cP9'][i][0],'o', color='orange')
    d.append(pro_data['cT'][i][0])
    t.append(i-24*5)

d = array(d)
t = array(t)
t = concatenate((t[range(41,72)],t[range(0,40)]))
d = concatenate((d[range(41,72)],d[range(0,40)]))    

plot(t,d,'-', color='orange', lw=3)   
#plot_dict(w_data,'cP9',sf=10000)
yticks(fontsize=16)
xlim(-24,24*2)
ylim(1e2,1e5)
xticks(range(0,24*2,12),fontsize=16)
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(24+12,24*2, color='gray',alpha=0.3)
axvspan(12,24, color='gray',alpha=0.3)
xlabel('Time in constant light (h)', fontsize=22)
ylabel('#monomers/cell', fontsize=22)
yscale('log')
tight_layout()
savefig('U2019_4_cT_rescaled_total.svg',format='svg',dpi=600, transparent=True)
savefig('U2019_4_cT_rescaled_total.pdf',format='pdf',dpi=600, transparent=True)
In [125]:
## print log(2)/U2017.m15
print log(2)/U2017.m23
print log(2)/U2017.m17
print log(2)/U2017.m24

0.4681601576364287
1.1426575868077768
9.040130030680263
In [208]:
temp_pro_data = []
temp_pro_data_time = []
for i in sorted(pro_data['cT'].keys()):
    temp_pro_data.append(pro_data['cT'][i][0])
    temp_pro_data_time.append(i)
In [209]:
#plot(temp_pro_data_time[:20],temp_pro_data[:20]/max(temp_pro_data[:20]),'o')
plot(temp_pro_data_time,temp_pro_data,'o')
print min(temp_pro_data[20:])
print max(temp_pro_data[20:])

452.94743129192386
29981.58067324995
In [151]:
figure(figsize=(10,6))
yticks(fontsize=16)
xticks(range(0,24*2,12),fontsize=16)
plot(res['time']-24*11,(res['[cP7]'])/max(res['[cP7]']),'k', label='cP7')
plot(res['time']-24*11,(res['[cL_m]']*150)/max(res['[cL_m]']*150),'b', label='cL_m')
rep = U2017.gp7**2*U2017.g1**2/(U2017.gp7**2*U2017.g1**2+res['[cP7]']**2)
plot(res['time']-24*11,1-rep, 'r', label='cP7 repression')
plot(array(temp_pro_data_time)-24*5,temp_pro_data/max(temp_pro_data),'o', color='orange', label='PRR7 SPD')
#plot_dict(w_data,'cP7',sf=50000)
title('cP7', fontsize=22)
xlim(-24,24*2)
ylim(0,1.7)
axvspan(-12,0, color='gray',alpha=0.7)
ylabel('Normalised to max', fontsize=22)
xlabel('Time in constant light (h)', fontsize=22)
axvspan(24+12,24*2, color='gray',alpha=0.1)
axvspan(12,24, color='gray',alpha=0.1)
legend(loc='upper right', fontsize=16)

temp_time =[]
temp_rna =[]
for i in sorted(network_data['col_0_LL']['log_cC_m'].keys()):
    temp_time.append(i)
    temp_rna.append(exp(network_data['col_0_LL']['log_cC_m'][i][0]))

plot(array(temp_time[:])-24*10,temp_rna[:]/max(temp_rna[:]),'bo-')

savefig('cL_rise.png', format='png',dpi=300)
In [152]:
temp_time =[]
temp_rna =[]
for i in sorted(network_data['col_0_LL']['log_cC_m'].keys()):
    temp_time.append(i)
    temp_rna.append(exp(network_data['col_0_LL']['log_cC_m'][i][0]))

plot(temp_time[1:],temp_rna[1:],'bo-')
Out[152]:
[<matplotlib.lines.Line2D at 0x12bb3a610>]
In [153]:
TOPLESS_expression = array([337.836,266.11,225.958,399.275,470.614,438.132,
                            320.583,342.415,279.855,428.996,411.956,324.224])
plot(range(0,24*2,4),TOPLESS_expression/max(TOPLESS_expression),'bo-')
Out[153]:
[<matplotlib.lines.Line2D at 0x129fd1c50>]
In [154]:
figure(figsize=(10,6))
yticks(fontsize=16)
xticks(range(0,24*2,12),fontsize=16)
plot(res['time']-24*11,(res['[cP7]'])/max(res['[cP7]']),'k', label='cP7')
plot(res['time']-24*11,(res['[cL_m]']*150)/max(res['[cL_m]']*150),'b', label='cL_m')
plot(range(0,24*2,4),TOPLESS_expression/max(TOPLESS_expression),'go-', label='TPL Diurnal Database')
plot(array(temp_pro_data_time)-24*5,temp_pro_data/max(temp_pro_data),'o', color='orange', label='PRR7 SPD')
#plot_dict(w_data,'cP7',sf=50000)
title('cP7', fontsize=22)
xlim(-24,24*2)
ylim(0,1.7)
axvspan(-12,0, color='gray',alpha=0.7)
ylabel('Normalised to max', fontsize=22)
xlabel('Time in constant light (h)', fontsize=22)
axvspan(24+12,24*2, color='gray',alpha=0.1)
axvspan(12,24, color='gray',alpha=0.1)
legend(loc='upper right', fontsize=16)
savefig('topless_progile.png', format='png',dpi=300)
In [155]:
U2017sf
Out[155]:
KeyedList([('gml', 9696.3284783153849), ('gm9', 3645.1361509496446), ('gm7', 1854.8827904077939), ('gm5', 530.27711530965303), ('gmt', 4881.4791683880294), ('gm4', 4919.2311439361465), ('gm3', 2847.6321058166104), ('gmg', 4987.7906816742825)])
In [211]:
U2017k = te.loada('Models/U2017kinase.txt')

U2017p = Utility.load('Models/params_4.bp')
U2017sf = Utility.load('Models/sfU2017.bp')
for i in U2017p.keys():
    U2017k.setValue(i,U2017p.getByKey(i))
for i in U2017sf.keys():
    U2017k.setValue(i,U2017sf.getByKey(i))
for i in p.keys():
    U2017k.setValue(i,p[i])
U2017k.setValue('k_strength',0)
reskbefore = U2017k.simulate(0,24*15,3000)
print log(2)/U2017k.m15
print log(2)/U2017k.m23
U2017k.setValue('k_strength',0.15)
U2017k.setValue('m15',0.3)
U2017k.setValue('m23',0.01)
print log(2)/0.3
print log(2)/0.01
resk = U2017k.simulate(0,24*15,3000)

1.004126352466796
0.39976022536822303
2.3104906018664844
69.31471805599453
In [216]:
figure(figsize=(10,6))
yticks(fontsize=16)
xticks(range(0,24*2,12),fontsize=16)
plot(res['time']-24*11,(resk['[cP7]']),'k', label='cP7')
plot(res['time']-24*11,(resk['[cL_m]']*150),'b', label='cL_m*150')
plot(reskbefore['time']-24*11,(reskbefore['[cL_m]']*150),'b--', label='cL_m*150 no K')
plot(res['time']-24*11,(resk['[cK]']*150),'r', label='HKinase')
#plot(range(0,24*2,4),TOPLESS_expression,'go-', label='TPL Diurnal Database')
#plot(array(temp_pro_data_time)-24*5,temp_pro_data,'o', color='orange', label='PRR7 SPD')
#plot_dict(w_data,'cP7',sf=50000)
#title('cP7', fontsize=22)
xlim(-24,24*2)
#ylim(0,1.7)
axvspan(-12,0, color='gray',alpha=0.7)
ylabel('#monomers/cell', fontsize=22)
xlabel('Time in constant light (h)', fontsize=22)
axvspan(24+12,24*2, color='gray',alpha=0.3)
axvspan(12,24, color='gray',alpha=0.3)
axvspan(12+24*2,24*3, color='gray',alpha=0.3)
temp_time =[]
temp_rna =[]
for i in sorted(network_data['col_0_LL']['log_cC_m'].keys()):
    temp_time.append(i)
    temp_rna.append(exp(network_data['col_0_LL']['log_cC_m'][i][0]))

plot(array(temp_time[:])-24*11,array(temp_rna[:])*1.3e2,'bo-', label='CCA1 mRNA TiMet')
legend(loc='upper right', fontsize=16)
yscale('log')
#savefig('hkinase_profile.png', format='png',dpi=300)
In [126]:
figure(dpi=600)
d = []
t = []
for i in pro_data['cE3'].keys():
    #plot(i-24*5,pro_data['cP9'][i][0],'o', color='orange')
    d.append(pro_data['cE3'][i][0])
    t.append(i-24*5)

d = array(d)
t = array(t)
t = concatenate((t[range(41,72)],t[range(0,40)]))
d = concatenate((d[range(41,72)],d[range(0,40)]))   
plot(t,d,'-', color='orange', lw=3) 

xticks(range(0,24*2,12),fontsize=16)
plot(res['time']-24*11,(res['[cE3c]']), 'k')
#for i in pro_data['cE3'].keys():
#    plot(i-24*5,pro_data['cE3'][i][0],'o', color='orange')
#plot_dict(w_data,'cE3',sf=30000)
yticks(fontsize=16)
xlim(-24,24*2)
ylim(1e3,1e5)
yscale('log')
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(24+12,24*2, color='gray',alpha=0.3)
axvspan(12,24, color='gray',alpha=0.3)
ylabel('#monomers/cell', fontsize=22)
xlabel('Time in constant light (h)', fontsize=22)
tight_layout()
savefig('U2019_4_cE3_rescaled.svg',format='svg',dpi=600,transparent=True)
savefig('U2019_4_cE3_rescaled.pdf',format='pdf',dpi=600,transparent=True)
In [127]:
figure(dpi=600)
d = []
t = []
for i in pro_data['cE4'].keys():
    #plot(i-24*5,pro_data['cP9'][i][0],'o', color='orange')
    d.append(pro_data['cE4'][i][0])
    t.append(i-24*5)

d = array(d)
t = array(t)
t = concatenate((t[range(41,72)],t[range(0,40)]))
d = concatenate((d[range(41,72)],d[range(0,40)]))   
plot(t,d,'-', color='orange', lw=3) 

xticks(range(0,24*2,12),fontsize=16)
plot(res['time']-24*11,(res['[cE4]']), 'k')
#for i in pro_data['cE3'].keys():
#    plot(i-24*5,pro_data['cE3'][i][0],'o', color='orange')
#plot_dict(w_data,'cE3',sf=30000)
yticks(fontsize=16)
xlim(-24,24*2)
ylim(1e2,1e6)
yscale('log')
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(24+12,24*2, color='gray',alpha=0.3)
axvspan(12,24, color='gray',alpha=0.3)
ylabel('#monomers/cell', fontsize=22)
xlabel('Time in constant light (h)', fontsize=22)
tight_layout()
savefig('U2019_4_cE4_rescaled.svg',format='svg',dpi=600,transparent=True)
savefig('U2019_4_cE4_rescaled.pdf',format='pdf',dpi=600,transparent=True)
In [128]:
figure(dpi=600)
d = []
t = []
for i in pro_data['cLUX'].keys():
    #plot(i-24*5,pro_data['cP9'][i][0],'o', color='orange')
    d.append(pro_data['cLUX'][i][0])
    t.append(i-24*5)

d = array(d)
t = array(t)
t = concatenate((t[range(41,72)],t[range(0,40)]))
d = concatenate((d[range(41,72)],d[range(0,40)]))   
plot(t,d,'-', color='orange', lw=3) 

xticks(range(0,24*2,12),fontsize=16)
plot(res['time']-24*11,(res['[cLUX]']), 'k')
#for i in pro_data['cE3'].keys():
#    plot(i-24*5,pro_data['cE3'][i][0],'o', color='orange')
#plot_dict(w_data,'cE3',sf=30000)
yticks(fontsize=16)
xlim(-24,24*2)
ylim(1e4,1e6)
yscale('log')
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(24+12,24*2, color='gray',alpha=0.3)
axvspan(12,24, color='gray',alpha=0.3)
ylabel('#monomers/cell', fontsize=22)
xlabel('Time in constant light (h)', fontsize=22)
tight_layout()
savefig('U2019_4_cLUX_rescaled.svg',format='svg',dpi=600,transparent=True)
savefig('U2019_4_cLUX_rescaled.pdf',format='pdf',dpi=600,transparent=True)
In [ ]:
U2

Translation rates

In [167]:
print max(get_max[20:])
print min(get_max[20:])

128663.174159
62252.9310291
In [240]:
pbm_affinities = read_excel('data/LUX_affinities_PBM_U2017.xlsx')
pbm_affinities
Out[240]:
Model Gene promoter Kd value
0 PBM LUX 0.260289 0.260289
1 PBM PRR9 0.276830 0.276830
2 PBM CCA1 0.393578 0.393578
3 PBM PRR5 1.977416 1.977416
4 PBM TOC1 0.685716 0.685716
5 U2017 PRR9 9.250000 0.010300
6 U2017 TOC1 10.200000 0.011500
7 U2017 LUX 8.850000 0.009900
8 U2017 GI 3.580000 0.004010
In [241]:
plt.style.use('seaborn-white')
figure(figsize=(10,6))
sns.set_palette("Paired")
pbm_vs_U2017 =  read_excel('data/LUX_affinities_PBM_U2017.xlsx')
sns.barplot(x='Gene promoter', y='Kd',data=pbm_vs_U2017, hue='Model')
legend(loc='upper left', fontsize=16)
ylabel('$k_d$ nM', fontsize=22)
xlabel('Gene promoter',fontsize=22)
ylim(0,12)
#yscale('log')
xticks(fontsize=16)
yticks(fontsize=16)
savefig('EC_Kd_comparision_absolute.png', format='png', dpi=300)


NameErrorTraceback (most recent call last)
<ipython-input-241-38351d91821a> in <module>()
      1 plt.style.use('seaborn-white')
      2 figure(figsize=(10,6))
----> 3 sns.set_palette("Paired")
      4 pbm_vs_U2017 =  read_excel('data/LUX_affinities_PBM_U2017.xlsx')
      5 sns.barplot(x='Gene promoter', y='Kd',data=pbm_vs_U2017, hue='Model')

NameError: name 'sns' is not defined
<Figure size 720x432 with 0 Axes>
In [242]:
U2017k = te.loada('Models/U2017KBOA.txt')

U2017p = Utility.load('Models/params_4.bp')
U2017sf = Utility.load('Models/sfU2017.bp')
for i in U2017p.keys():
    U2017k.setValue(i,U2017p.getByKey(i))
for i in U2017sf.keys():
    U2017k.setValue(i,U2017sf.getByKey(i))
for i in p.keys():
    U2017k.setValue(i,p[i])


U2017k.setValue('k_strength',0.1)
U2017k.setValue('pBOA',0.0)
U2017k.setValue('pBOA',20)

U2017k.setValue('pBOA2',0.5)
U2017k.setValue('pBOA2',200)
#U2017k.setValue('p27',0)
U2017k.setValue('m15',0.3)
U2017k.setValue('m23',0.3)
#U2017k.setValue('aff',0.1)
U2017k.setValue('aff',0.01)
U2017k.setValue('fass_ev',200)
#U2017k.setValue('fass_ev',1)
print log(2)/0.3
print log(2)/0.01
resk = U2017k.simulate(0,24*15,3000)

figure(figsize=(10,6))
yticks(fontsize=16)
xticks(range(0,24*4,12),fontsize=16)
plot(res['time']-24*11,(resk['[cBOA]']),'g', label='cBOA')
#plot(res['time']-24*11,(resk['[cEC]']),'r-', label='cEC')
plot(res['time']-24*11,(resk['[cL_m]']),'b-', label='BOA-OX cL_m')
xlim(-24,24*10)
#ylim(0,1.7)
axvspan(-12,0, color='gray',alpha=0.7)
ylabel('log(#molecules/cell)', fontsize=22)
xlabel('Time in constant light (h)', fontsize=22)
axvspan(24+12,24*2, color='gray',alpha=0.1)
axvspan(12,24, color='gray',alpha=0.1)
axvspan(12+24*2,24*3, color='gray',alpha=0.1)
yscale('log')


title('B', fontsize=22)
legend(loc='upper right', fontsize=16)
#savefig('BOA_intro.png',format='png',dpi=300)

2.3104906018664844
69.31471805599453
Out[242]:
<matplotlib.legend.Legend at 0x7f5e359a8910>
In [274]:
amin_sec = 3.
rib_tranc = 10.
PRR7_length = 468.
PRR7_trans_rate = 1/(PRR7_length/amin_sec/10)*3600
plot(res['time']-24*11,res['[cP9_m]']*PRR7_trans_rate,'k', label='Theoretical Max trans')
plot(res['time']-24*11,U2017.getValue('gp7')*U2017.getValue('p8')*res['[cP9_m]']/U2017.getValue('gm9'), 'b--',label='Rescaled')
xticks(range(0,24*2,12), fontsize=22)
yticks(fontsize=22)
xlim(-24,24*2)
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(24+12,24*2, color='gray',alpha=0.3)
axvspan(12,24, color='gray',alpha=0.3)
legend(loc='lower right', fontsize=16)
ylabel('#molecules/h', fontsize=22)
yscale('log')
In [270]:
figure(figsize=(12,10), dpi=600)

subplot(221)
amin_sec = 3.
rib_tranc = 10.
PRR7_length = 468.
PRR7_trans_rate = 1/(PRR7_length/amin_sec/10)*3600
plot(res['time']-24*11,res['[cP9_m]']*PRR7_trans_rate,'k', label='Theoretical Max trans')
plot(res['time']-24*11,U2017.getValue('gp7')*U2017.getValue('p8')*res['[cP9_m]']/U2017.getValue('gm9'), 'b--',label='Rescaled')
yscale('log')
xticks(range(0,24*2,12), fontsize=22)
yticks(fontsize=22)
xlim(-24,24*2)
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(24+12,24*2, color='gray',alpha=0.3)
axvspan(12,24, color='gray',alpha=0.3)
legend(loc='lower right', fontsize=16)
ylabel('#molecules/h', fontsize=22)

title('cP9', fontsize=22)

subplot(222)
amin_sec = 3.
rib_tranc = 10.
PRR7_length = 618.
PRR7_trans_rate = 1/(PRR7_length/amin_sec/10)*3600
plot(res['time']-24*11,res['[cP7_m]']*PRR7_trans_rate,'k', label='Theoretical max trans.')
plot(res['time']-24*11,U2017.getValue('gp7')*U2017.getValue('p9')*res['[cP7_m]']/U2017.getValue('gm7'), 'b--',label='Rescaled')
yscale('log')
xticks(range(0,24*2,12), fontsize=22)
yticks(fontsize=22)
xlim(-24,24*2)
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(24+12,24*2, color='gray',alpha=0.3)
axvspan(12,24, color='gray',alpha=0.3)


title('cP7', fontsize=22)

subplot(223)
amin_sec = 3.
rib_tranc = 10.
PRR7_length = 558.
PRR7_trans_rate = 1/(PRR7_length/amin_sec/10)*3600
plot(res['time']-24*11,res['[cP5_m]']*PRR7_trans_rate,'k', label='Theoretical Max trans')
plot(res['time']-24*11,U2017.getValue('gp7')*U2017.getValue('p10')*res['[cP5_m]']/U2017.getValue('gm5'),'b--', label='U2017 trans')
yscale('log')
xticks(range(0,24*2,12), fontsize=22)
yticks(fontsize=22)
xlim(-24,24*2)
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(24+12,24*2, color='gray',alpha=0.3)
axvspan(12,24, color='gray',alpha=0.3)

ylabel('#molecules/h', fontsize=22)
xlabel('Time in constant light (h)', fontsize=22)
title('cP5', fontsize=22)

subplot(224)
amin_sec = 3.
rib_tranc = 10.
PRR7_length = 618.
PRR7_trans_rate = 1/(PRR7_length/amin_sec/10)*3600
plot(res['time']-24*11,res['[cT_m]']*PRR7_trans_rate,'k', label='Theoretical Max trans')
plot(res['time']-24*11,U2017.getValue('gp7')*U2017.getValue('p4')*res['[cT_m]']/U2017.getValue('gmt'), 'b--',label='Rescaled')
yscale('log')
xticks(range(0,24*2,12), fontsize=22)
yticks(fontsize=22)
xlim(-24,24*2)
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(24+12,24*2, color='gray',alpha=0.3)
axvspan(12,24, color='gray',alpha=0.3)


xlabel('Time in constant light (h)', fontsize=22)
title('cT', fontsize=22)
tight_layout()
savefig('PRRS__translation_rates.svg', format='svg', dpi=600, transparent=True)
In [259]:
from scipy.constants import Avogadro
In [260]:
print U2017.getValue('g9')
mean(array([U2017.getValue('g9'),
            U2017.getValue('g10'),
            U2017.getValue('g5'),
            U2017.getValue('g6'),
            U2017.getValue('g16'),
            U2017.getValue('g15'),
            U2017.getValue('g12')]))

0.267310750786
Out[260]:
0.29790048402767483
In [261]:
print 'g9\t', U2017.getValue('gl')*(U2017.getValue('p1')+U2017.getValue('p2'))/U2017.getValue('gml')

g9	158.249727823
In [262]:
U2017.getValue('m3')*2
Out[262]:
0.40110991504829296
In [263]:
A = 0.1168e-12
B = 0.112533e-12
(U2017.getValue('g9')*U2017.getValue('gl')/1100)/(Avogadro*.112533e-12)
Out[263]:
3.0448506487979457e-10
In [264]:
print 'g9\t', U2017.getValue('g9')*U2017.getValue('gl')/(Avogadro*.112533e-12)*1e9
print 'g10\t', U2017.getValue('g10')*U2017.getValue('gl')/(Avogadro*0.112533e-12)*1e9
print 'g12\t', U2017.getValue('g12')*U2017.getValue('gl')/(Avogadro*0.112533e-12)*1e9
print 'g5\t', U2017.getValue('g5')*U2017.getValue('gl')/(Avogadro*0.112533e-12)*1e9
print 'g6\t', U2017.getValue('g6')*U2017.getValue('gl')/(Avogadro*0.112533e-12)*1e9
print 'g16\t', U2017.getValue('g16')*U2017.getValue('gl')/(Avogadro*0.112533e-12)*1e9
print 'g15\t', U2017.getValue('g15')*U2017.getValue('gl')/(Avogadro*0.112533e-12)*1e9

g9	334.933571368
g10	619.642572909
g12	271.418814406
g5	179.422681504
g6	353.102873517
g16	383.312968504
g15	470.998621027
In [265]:
PBM_affinties = array([0.34,
0.25,
0.57,
0.44,
0.70,
0.40,
0.31,
0.40])



print PBM_affinties

PBM_affinties = array([0.44,
0.31,
0.47,
0.4,
0.25,
0.34])

U2017_cca1_aff  = array([U2017.getValue('g9')*U2017.getValue('gl')/(Avogadro*A)*1e9/1100,
U2017.getValue('g10')*U2017.getValue('gl')/(Avogadro*A)*1e9/1100,
U2017.getValue('g12')*U2017.getValue('gl')/(Avogadro*A)*1e9/1100,
U2017.getValue('g5')*U2017.getValue('gl')/(Avogadro*A)*1e9/1100,
U2017.getValue('g6')*U2017.getValue('gl')/(Avogadro*A)*1e9/1100,
U2017.getValue('g16')*U2017.getValue('gl')/(Avogadro*A)*1e9/1100,
U2017.getValue('g15')*U2017.getValue('gl')/(Avogadro*A)*1e9/1000])

U2017_cca1_aff  = U2017_cca1_aff
print U2017_cca1_aff*1100
print U2017_cca1_aff
print PBM_affinties

[0.34 0.25 0.57 0.44 0.7  0.4  0.31 0.4 ]
[322.6975992  597.00545939 261.50319727 172.86791625 340.2031307
 369.30957436 499.17103255]
[0.29336145 0.54273224 0.23773018 0.15715265 0.30927557 0.33573598
 0.45379185]
[0.44 0.31 0.47 0.4  0.25 0.34]
In [266]:
import seaborn as sns


ImportErrorTraceback (most recent call last)
<ipython-input-266-a84c0541e888> in <module>()
----> 1 import seaborn as sns

ImportError: No module named seaborn
In [267]:
pbm_vs_U2017 =  read_excel('data/PBM_vs_U20172.xlsx')
pbm_vs_U2017
Out[267]:
Gene promoter Model value ratio
0 PRR9 U2017.2 0.440000 0.909091
1 PRR7 U2017.2 0.310000 1.290323
2 PRR5 U2017.2 0.470000 0.851064
3 TOC1 U2017.2 0.400000 1.000000
4 LUX U2017.2 0.340000 1.176471
5 ELF4 U2017.2 0.250000 1.600000
6 GI U2017.2 0.570000 0.701754
7 PRR9 EMA selected 0.293360 0.535696
8 PRR7 EMA selected 0.542730 0.289558
9 PRR5 EMA selected 0.237729 0.661055
10 TOC1 EMA selected 0.157152 1.000000
11 LUX EMA selected 0.309274 0.508131
12 ELF4 EMA selected 0.309274 0.508131
13 GI EMA selected 0.453790 0.346310
In [183]:
plt.style.use('seaborn-white')
figure(figsize=(10,6))
sns.set_palette("Paired")
pbm_vs_U2017 =  read_excel('data/PBM_vs_U20172.xlsx')
sns.barplot(x='Gene promoter', y='ratio',data=pbm_vs_U2017, hue='Model')
legend(loc='upper left', fontsize=16)
ylabel('$k_d/k_d^{TOC1}$', fontsize=22)
xlabel('Gene promoter',fontsize=22)
title('B', fontsize=22)
xticks(fontsize=16)
yticks(fontsize=16)
savefig('Kd_comparision.pdf', format='pdf', dpi=300)
In [280]:
U2017 = te.loada('Models/U2017.txt')
#
U2017.reset() #= te.loada('Models/U2017.txt')
for i in U2017p.keys():
    U2017.setValue(i,U2017p.getByKey(i))

U2017.setValue('ge3',1)
range(0,24*10,1)
res = U2017.simulate(0,24*15,24*15)
#plot(res['time'],res['[cLUX]']/max(res['[cLUX]']))
U2017 = te.loada('Models/U2017.txt')

U2017p = Utility.load('Models/params_4.bp')
U2017sf = Utility.load('Models/sfU2017.bp')
for i in U2017p.keys():
    U2017.setValue(i,U2017p.getByKey(i))
for i in U2017sf.keys():
    U2017.setValue(i,U2017sf.getByKey(i))
for i in p.keys():
    U2017.setValue(i,p[i])

res = U2017.simulate(0,24*15,3000)


figure(figsize=(15,5))

subplot(121)
plot(res['time']-24*11,res['[cG_m]'],'k--', label='cG_m')
print 'G', max(res['[cG_m]'][1500:])
#yscale('log')
xticks(range(0,24*2,12), fontsize=16)
yticks(fontsize=16)
xlim(-24,24*2)
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(24+12,24*2, color='gray',alpha=0.3)
axvspan(12,24, color='gray',alpha=0.3)
legend(loc='upper right', fontsize=16)
ylabel('#Transcripts/cell', fontsize=22)
xlabel('Time in constant light (h)', fontsize=22)

ylim(0,400)
subplot(122)
plot(res['time']-24*11,res['[cT_m]'],'b--',label='cT_m')
print 'T', max(res['[cT_m]'][1500:])



U2017.setValue('g15',0.57)
U2017.setValue('g5',0.4)
res = U2017.simulate(0,24*15,3000)

subplot(121)
plot(res['time']-24*11,res['[cG_m]'],'k', label='cG_m')
print 'GI', max(res['[cG_m]'][1500:])
legend(loc='upper right', fontsize=16)
subplot(122)
plot(res['time']-24*11,res['[cT_m]'],'b',label='cT_m')
legend(loc='upper right', fontsize=16)
print 'TOC', max(res['[cT_m]'][1500:])


#yscale('log')
xticks(range(0,24*2,12), fontsize=16)
yticks(fontsize=16)
ylim(0,400)
xlim(-24,24*2)
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(24+12,24*2, color='gray',alpha=0.3)
axvspan(12,24, color='gray',alpha=0.3)
ylabel('#Transcripts/cell', fontsize=22)
xlabel('Time in constant light (h)', fontsize=22)

#savefig('GI_and_TOC1_params_updated.png', format='png',dpi=300)

G 237.70860756406006
T 71.54658805989442
GI 292.44079777486814
TOC 132.0887656011767
Out[280]:
Text(0.5,0,'Time in constant light (h)')
In [ ]:
print 'g8\t', U2017.getValue('g8')*U2017.getValue('ge3')/(Avogadro*A)*1e9
print 'g4\t', U2017.getValue('g4')*U2017.getValue('ge3')/(Avogadro*A)*1e9
print 'g2\t', U2017.getValue('g2')*U2017.getValue('ge3')/(Avogadro*A)*1e9
print 'g14\t', U2017.getValue('g14')*U2017.getValue('ge3')/(Avogadro*A)*1e9
In [ ]:
print 'g1\t', U2017.getValue('g1')*U2017.getValue('gp7')/(Avogadro*A)*1e9
print 'g7\t', U2017.getValue('g7')*U2017.getValue('gp7')/(Avogadro*A)*1e9
print 'g13\t', U2017.getValue('g13')*U2017.getValue('gp7')/(Avogadro*A)*1e9
In [ ]:
print 'cP9m\tg8\t',U2017.getValue('ge3')*U2017.getValue('g8')/(Avogadro*A)
print 'cTm\tg4\t',U2017.getValue('ge3')*U2017.getValue('g4')/(Avogadro*A)
print 'cE4m\tg2\t',U2017.getValue('ge3')*U2017.getValue('g2')/(Avogadro*A)
print 'cGm\tg14\t',U2017.getValue('ge3')*U2017.getValue('g14')/(Avogadro*A)

Evening complex

In [ ]:
mean(array([0.195076690644,
0.360900885544,
0.158083538406,
0.104501865224,
0.20565910948,
0.223254495126,
0.301758157905]))

PBM CCA1 concentration

In [210]:
figure(figsize=(7,5))
x = linspace(0,20,100)
used = 1e-6/(42.5+66.975)*1e9
plot(x,x**2/(x**2+2.5**2), label = 'CCR2')
plot(x,x**2/(x**2+5**2), label='wt PRR9')
plot(x,x**2/(x**2+11**2), label='PRR9 CBS')
axvline(used, color='red', label='PBM concentration')
ylabel('P(CCA1 bound)', fontsize=22)
xlabel('MPB-CCA1 nM', fontsize=22)
yticks(fontsize=22)
xticks(fontsize=22)
tight_layout()
legend(loc='lower right', fontsize=16)
savefig('CCA1_PBM_affinity.png', format='png', dpi=300)
In [ ]:
figure(figsize=(15,5))
subplot(121)
AT = linspace(1e-9,1e-6,1000000)
BT = 1100/(Avogadro*1.12e-13)
k=0.35*1e-9
A = 0.5*(AT-BT-k+sqrt((AT-BT-k)**2+4*AT*k))
plot(AT,A,'blue')
plot(AT,AT,'red')
yscale('log')
xscale('log')
#ylim(1e-2,1e5)
#xlim(1e0,1e5)
xticks(fontsize=22)
yticks(fontsize=22)
ylabel('Active cL nM', fontsize=22)
xlabel('nM', fontsize=22)
subplot(122)
prr9_kd = 0.355309957403*1e-9
plot(AT,A/(prr9_kd+A), color='blue', label='Titration')
plot(AT,AT**2/(prr9_kd**2+AT**2), color='red', label='Graded')
xscale('log')
xlabel('nM')
#xlim(1,1e3)
ylim(0,1.2)
legend(loc='lower right', fontsize=16)
xticks(fontsize=22)
yticks(fontsize=22)
ylabel('cL binding probability', fontsize=22)
xlabel('nM', fontsize=22)
In [ ]:
cL = res['[cL]']/(Avogadro*1.12e-13)
plot(res['time'],cL)
In [ ]:
(U2017.getValue('gl')*U2017.getValue('g9'))/(Avogadro*A)
In [ ]:
A = 0.1168e-12
#plot(res['time'],res['[cL]']**2/(k**2+res['[cL]']**2),'o')
cL_con = 0.5*(cL-BT-k+sqrt((cL-BT-k)**2+4*cL*k))
plot(res['time'],cL_con/((U2017.getValue('gl')*U2017.getValue('g9')/(Avogadro*A)/1100+cL_con)),'-')
#yscale('log')
In [ ]:
U2017.getValue('gl')*U2017.getValue('g9')/(Avogadro*A)
In [ ]:
figure(figsize=(10,6))
xticks(range(0,24*2,12),fontsize=16)
plot(res['time']-24*11,res['[cEC]'], 'k')

yticks(fontsize=16)
xlim(-24,24*2)
#ylim(0,3e4)
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(24+12,24*2, color='gray',alpha=0.1)
axvspan(12,24, color='gray',alpha=0.1)
ylabel('#monomers/cell', fontsize=22)
xlabel('Time in constant light (h)', fontsize=22)
title('cEC', fontsize=22)
savefig('EC_dynamics_rescaled.png', format='png',dpi=300)
In [ ]:
figure(figsize=(10,6))
xticks(range(0,24*2,12),fontsize=16)
plot(res['time']-24*11,res['[cT_m]']*1e2, 'k')
plot(res['time']-24*11,res['[cEC]'], 'k')



yticks(fontsize=16)
xlim(-24,24*2)
#ylim(0,3e4)
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(24+12,24*2, color='gray',alpha=0.1)
axvspan(12,24, color='gray',alpha=0.1)
ylabel('#monomers/cell', fontsize=22)
xlabel('Time in constant light (h)', fontsize=22)
title('cEC', fontsize=22)
In [ ]:
figure(figsize=(10,6))
xticks(range(0,24*2,12),fontsize=16)
plot(res['time']-24*11,res['[cT_m]']*1e2, 'k')
plot(res['time']-24*11,res['[cEC]'], 'k')



yticks(fontsize=16)
xlim(-24,24*2)
#ylim(0,3e4)
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(24+12,24*2, color='gray',alpha=0.1)
axvspan(12,24, color='gray',alpha=0.1)
ylabel('#monomers/cell', fontsize=22)
xlabel('Time in constant light (h)', fontsize=22)
title('cEC', fontsize=22)
In [ ]:
figure(figsize=(10,6))
xticks(range(0,24*2,12),fontsize=16)
U2017.reset()
res = U2017.simulate(0,24*15,3000)
plot(res['time']-24*11,res['[cT_m]']*1e2, 'k')
plot(res['time']-24*11,res['[cEC]'], 'k')



yticks(fontsize=16)
xlim(-24,24*2)
#ylim(0,3e4)
axvspan(-12,0, color='gray',alpha=0.7)
axvspan(24+12,24*2, color='gray',alpha=0.1)
axvspan(12,24, color='gray',alpha=0.1)
ylabel('#monomers/cell', fontsize=22)
xlabel('Time in constant light (h)', fontsize=22)
title('cEC', fontsize=22)
In [ ]:
U2017 = te.loada('Models/U2017.txt')
#

U2017 = te.loada('Models/U2017.txt')

U2017p = Utility.load('Models/params_4.bp')
U2017sf = Utility.load('Models/sfU2017.bp')
for i in U2017p.keys():
    U2017.setValue(i,U2017p.getByKey(i))
for i in U2017sf.keys():
    U2017.setValue(i,U2017sf.getByKey(i))
for i in p.keys():
    U2017.setValue(i,p[i])

res = U2017.simulate(0,24*15,3000)
In [ ]:
cL_temp = linspace(1,1e3)
hill = U2017.p3*U2017.gl*cL_temp**2/(U2017.gl*U2017.g3**2+cL_temp**2)+cL_temp*U2017.m3
linear = [res['[cL_m]'][argmax(res['[cL]'])]*(U2017.p1+U2017.p2)/U2017.gml]*len(cL_temp)
In [ ]:
plot(cL_temp,hill)
plot(cL_temp,linear)
xscale('log')
In [ ]:
k = []
for i in pro_data['cL'].keys():
    k.append(pro_data['cL'][i][0])
k=array(k)
In [ ]:
#in aminoacids 
amin_sec = 3.
rib_tranc = 10.
CCA1_length = 609.
CCA1_trans_rate = 1/(CCA1_length/amin_sec/10)*3600
print CCA1_trans_rate
amin_sec = 3.
rib_tranc = 10.
LHY_length = 646.
LHY_trans_rate = 1/(LHY_length/amin_sec/10)*3600
print LHY_trans_rate*exp(lhyinter(sorted(pro_data['cL'].keys())))
In [1208]:
figure(figsize=(7,5))
p = []
for i in pro_data['cL'].keys():
    p.append(pro_data['cL'][i][0])

    
synth =LHY_trans_rate*exp(lhyinter(sorted(pro_data['cL'].keys())))

deg = (1.76)*array(p)

plot(array(sorted(pro_data['cL'].keys()))-24*4-24,(synth/deg), 'b',label='SM')
xlim(0,24*3)

synth = res['[cL_m]']*(U2017.p1+U2017.p2)*U2017.gl/U2017.gml
deg = U2017.m2*res['[cL]']+U2017.p3*U2017.gl*res['[cL]']**2/((U2017.gl+U2017.g3)**2+res['[cL]']**2)

plot(res['time']-24*10-24,(synth/deg), 'g',label='U2017')
xlim(0,24*4)
xticks(range(0,24*3,12),fontsize=16)
yticks(fontsize=16)
xlim(-24,24*2)
yscale('log')
axvspan(12-24,24-24,0, color='gray',alpha=0.7)
axvspan(12+24-24,24*2-24, color='gray',alpha=0.1)
axvspan(12+24*2-24,24*3-24, color='gray',alpha=0.1)
ylabel('log(synt/deg)', fontsize=22)
xlabel('Time in constant light (h)', fontsize=22)
#title('A', fontsize=22)
legend(loc='upper right', fontsize=16)
savefig('cL_synth_deg.png', format='png',dpi=300)
In [ ]:
print log10(max(res['[cL]'][1000:1300]))
print log10(min(res['[cL]'][1000:1300]))
In [ ]:
synth = res['[cL_m]']*(U2017.p1+U2017.p2)*U2017.gl/U2017.gml
deg = U2017.m2*res['[cL]']+U2017.p3*U2017.gl*res['[cL]']**2/((U2017.gl+U2017.g3)**2+res['[cL]']**2)
In [ ]:
plot(res['time']-24*10,(synth-deg))
xlim(0,24*4)
In [ ]:
In [ ]: