#Model fitting
library(minpack.lm)
library(bbmle)
library(plyr)
library(ggplot2)
#Observed data between 0uM and 10uM
x<-c(0,0,0,0.5,0.5,0.5,1,1,1,2.5,2.5,2.5,5,5,5,7.5,7.5,7.5,10,10,10)
y<-c(0.002857143,0.00152381,0.00212963,0.011598174,0.009054726,0.011746032,0.021904762,0.00854902,0.016464646,0.038600823,0.043166667,0.043508772,0.045142857,0.045367965,0.035365854,0.043474178,0.044122807,0.04971831,0.03822807,0.045063291,0.043288889)
#0.043508772<=0.053508772, 0.04971831<=0.05971831, 0.03822807<=0.033?
df<-data.frame(x, y)
#HILL FUNCTION
Hill<-nlsLM(y~((a+b*x^d)/(c+x^d)), start=list(a=0.002, b=0.044, c=1, d=2), trace=TRUE)
#LOGISYIC FUNCTION
Logt<-nlsLM(y~(a+b/(1+exp(4*c*(d-x)/b+2))), start=list(a=0.001, b=0.044, c=0.025, d=0.1), trace=TRUE)
#GOMPERTZ FUNCTION
Gomp<-nlsLM(y~(a+b*exp(-exp(c*2.71828*(d-x)/b+1))), start=list(a=0.001, b=0.044, c=0.025, d=0.1), trace=TRUE)
AICtab(Hill, Logt, Gomp)
# dAIC df
#Logt 0.0 5
#Gomp 1.8 5
#Hill 3.6 5
summary(Logt)
#Formula: y ~ (a + b/(1 + exp(4 * c * (d - x)/b + 2)))
#Parameters:
# Estimate Std. Error t value Pr(>|t|)
#a 0.0001397 0.0048826 0.029 0.977507
#b 0.0434417 0.0053800 8.075 3.22e-07 ***
#c 0.0222537 0.0054435 4.088 0.000766 ***
#d 0.2392522 0.3978514 0.601 0.555535
#Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#Residual standard error: 0.004108 on 17 degrees of freedom
#Number of iterations to convergence: 9
#Achieved convergence tolerance: 1.49e-08
plot(x, y, lwd=2, xlab=\'AuCl4_uM\', ylab=\'Relative_bGal\')
xfit <- seq(0, 10, 0.01)
yfit <- (0.0001397+0.041717/(1+exp(4*0.0222537*(0.2392522-xfit)/ 0.0434417+2)))
lines(spline(xfit, yfit), lwd=2)
#The effect of gold concentration
iGEM<-function(t, state, parms) {
with(as.list(c(state, parms)), {
if(Uin<0) {Uin<-0}
else {Uin<-Uin}
if(Uout<0) {Uout<-0}
else {Uout<-Uout}
dR <- (a + b/(1 + exp(4 * c * (d - Uin)/b + 2)))-R
dUin <- (e*(Uout-Uin))
dUout <- (e*(Uin-Uout))
list(c(dR, dUin, dUout))
})
}
library(deSolve)
parms<-c(
a=0.0001397,
b=0.0434417,
c=0.0222537,
d=0.2392522,
e=0.01)
times<-seq(0, 120, 0.05)
state<-c(
R=0.003464946,
Uin=0,
Uout=0)
out0 <- ode(y = state, times = times, func = iGEM, parms = parms)
state<-c(
R=0.003464946,
Uin=0,
Uout=0.5)
out1 <- ode(y = state, times = times, func = iGEM, parms = parms)
state<-c(
R=0.003464946,
Uin=0,
Uout=1)
out2 <- ode(y = state, times = times, func = iGEM, parms = parms)
state<-c(
R=0.003464946,
Uin=0,
Uout=2)
out3 <- ode(y = state, times = times, func = iGEM, parms = parms)
state<-c(
R=0.003464946,
Uin=0,
Uout=4)
out4 <- ode(y = state, times = times, func = iGEM, parms = parms)
state<-c(
R=0.003464946,
Uin=0,
Uout=6)
out5 <- ode(y = state, times = times, func = iGEM, parms = parms)
state<-c(
R=0.003464946,
Uin=0,
Uout=8)
out6 <- ode(y = state, times = times, func = iGEM, parms = parms)
state<-c(
R=0.003464946,
Uin=0,
Uout=10)
out7 <- ode(y = state, times = times, func = iGEM, parms = parms)
out<-rbind(out0, out1, out2, out3, out4, out5, out6, out7)
plot(out, pch=\'.\')
#The effect of membrane permeability
iGEM<-function(t, state, parms) {
with(as.list(c(state, parms)), {
if(Uin<0) {Uin<-0}
else {Uin<-Uin}
if(Uout<0) {Uout<-0}
else {Uout<-Uout}
dR <- (a + b/(1 + exp(4 * c * (d - Uin)/b + 2)))-R
dUin <- (e*(Uout-Uin))
dUout <- (e*(Uin-Uout))
list(c(dR, dUin, dUout))
})
}
time<-seq(0, 400, 1)
library(deSolve)
parms<-c(
a=0.0001397,
b=0.0434417,
c=0.0222537,
d=0.2392522,
e=0.01)
state<-c(
R=0.003464946,
Uin=0,
Uout=0)
out0 <- ode(y = state, times = times, func = iGEM, parms = parms)
times<-seq(0, 400, 0.1)
parms<-c(
a=0.0001397,
b=0.0434417,
c=0.0222537,
d=0.2392522,
e=0.01)
state<-c(
R=0.003464946,
Uin=0,
Uout=5)
out1 <- ode(y = state, times = times, func = iGEM, parms = parms)
parms<-c(
a=0.0001397,
b=0.0434417,
c=0.0222537,
d=0.2392522,
e=0.05)
state<-c(
R=0.003464946,
Uin=0,
Uout=5)
out2 <- ode(y = state, times = times, func = iGEM, parms = parms)
parms<-c(
a=0.0001397,
b=0.0434417,
c=0.0222537,
d=0.2392522,
e=0.005)
state<-c(
R=0.003464946,
Uin=0,
Uout=5)
out3 <- ode(y = state, times = times, func = iGEM, parms = parms)
out<-rbind(out0, out1, out2, out3)
plot(out, pch=\'.\')
#The effect of gold mineralization peptides
#mineralisastion to peptide degradation ratio
iGEM<-function(t, state, parms) {
with(as.list(c(state, parms)), {
if(Uin<0) {Uin<-0}
else {Uin<-Uin}
if(Uout<0) {Uout<-0}
else {Uout<-Uout}
if(Uout==0) {P<-0}
else {P<-P}
dR <- (a + b/(1 + exp(4 * c * (d - Uin)/b + 2)))-R
dUin <- (e*(Uout-Uin))
dUout <- (e*(Uin-Uout)-f*P)
dP<-(g*R/(h+R)-i*P)
list(c(dR, dUin, dUout, dP))
})
}
library(deSolve)
times<-seq(0, 150, 0.1)
parms<-c(
a=0.0001397,
b=0.0434417,
c=0.0222537,
d=0.2392522,
e=0.1,
f=0.1,
g=0.1,
h=0.1,
i=0.1)
state<-c(
R=0.003464946,
Uin=0,
Uout=0,
P=0)
out0 <- ode(y = state, times = times, func = iGEM, parms = parms)
state<-c(
R=0.003464946,
Uin=0,
Uout=5,
P=0)
out1 <- ode(y = state, times = times, func = iGEM, parms = parms)
parms<-c(
a=0.0001397,
b=0.0434417,
c=0.0222537,
d=0.2392522,
e=0.1,
f=0.12,
g=1,
h=1,
i=0.1)
out2 <- ode(y = state, times = times, func = iGEM, parms = parms)
parms<-c(
a=0.0001397,
b=0.0434417,
c=0.0222537,
d=0.2392522,
e=0.1,
f=0.2,
g=1,
h=1,
i=0.1)
out3 <- ode(y = state, times = times, func = iGEM, parms = parms)
out<-rbind(out0, out1, out2, out3)
plot(out, pch=\'.\')