dev/testingCode.R

In ffeng23/ADASPR: An R package to characterize the distributed affinities in a multi-ligand system

################taylor expansion
### to test taylor expansion of 1/(1-x) with polynomial of degree of 10.
x<-seq(0,1,by=0.01)
y<-1/(1-x[-1])
degree<-10
yexp<-rep(0,length(y))
for(i in c(1:length(y)))
{
	for(j in 1:degree)
	{
		yexp[i]<-yexp[i]+(x[i+1])^j
	}
}
yexp<-yexp+1
plot(x[-1],y)
lines(x[-1],yexp,col=2,lty=2)

######=------>
library(ADASPR)
#-------->testing code for multiligandmodel
kon<-c(2E3, 1E3, 3E3)
koff<-c(0.001, 0.01,0.05)

analyteConcentrations<-c(
							#5E-4, 4E-4, 3E-4, 2E-4, 
							1E-6, 2E-6, 1E-5)
associationLength<-2000 
dissociationLength<-2000
Rmax<-c(50,80,60)

mlgm<- new("MultiLigandModel", kon=kon, koff=koff, analyteConcentrations=analyteConcentrations, 
		associationLength=associationLength, dissociationLength=dissociationLength, Rmax=Rmax
		)
sData<-Simulate(mlgm)	
sData[[1]]
		
lgm1<-LangmuirModel(2e3, 0.001, analyteConcentrations=c(1e-6, 2e-6, 1e-5), Rmax=50, associationLength=100, dissociationLength=100)	
		
lgm2<-LangmuirModel(1e3, 0.01, analyteConcentrations=c(1e-6, 2e-6, 1e-5), Rmax=80, associationLength=100, dissociationLength=100)	

lgm3<-LangmuirModel(3e3, 0.05, analyteConcentrations=c(1e-6, 2e-6, 1e-5), Rmax=60, associationLength=100, dissociationLength=100)	

sData1<-Simulate(lgm1)
sData2<-Simulate(lgm2)
sData3<-Simulate(lgm3)
sDataSum<-GetObservedRUs(sData1[[1]],sData2[[1]])
sDataSum3<-GetObservedRUs(sDataSum,sData3[[1]])
#plotting R_AB
	plot(sData[[1]])
#plotting R_A
plot(sData[[2]])

library(ADASPR)
#------>for steady state spr
kon<-c(2E3, 1E3, 3E3)
koff<-c(0.001, 0.01,0.05)

analyteConcentrations<-c(
							#5E-4, 4E-4, 3E-4, 2E-4, 
							1E-4, 9E-5,6E-5, 3E-5, 2E-5)
associationLength<-1000
dissociationLength<-1000
Rmax<-c(50,80,60)

mlgm<- new("MultiLigandModel", kon=kon, koff=koff, analyteConcentrations=analyteConcentrations, 
		associationLength=associationLength, dissociationLength=dissociationLength, Rmax=Rmax
		)
sData<-Simulate(mlgm)		

#plotting R_AB
	plot(sData[[1]])
#plotting R_A
#plot(sData[[2]])	
e_k<-sum(koff/kon*(Rmax/sum(Rmax)))
e_k2<-sum((koff/kon)^2*(Rmax/sum(Rmax)))
e_k3<-sum((koff/kon)^3*(Rmax/sum(Rmax)))
e_k4<-sum((koff/kon)^4*(Rmax/sum(Rmax)))

FitSteadyStateSPR(sData[[1]], degree=5, steadyStateStart=950,steadyStateEnd=1000)

###-------->
###now testing 
#for non-steady state spr
kon<-c(2E3, 1E3, 3E3)
koff<-c(0.001, 0.01,0.05)

analyteConcentrations<-c(
							#5E-4, 4E-4, 3E-4, 2E-4, 
							1E-4, 9E-5,6E-5, 3E-5, 2E-5)
associationLength<-100 
dissociationLength<-100
Rmax<-c(50,80,60)

#mlgm<- new("MultiLigandModel", kon=kon, koff=koff, analyteConcentrations=analyteConcentrations, 
#		associationLength=associationLength, dissociationLength=dissociationLength, Rmax=Rmax
#		)
#sData<-Simulate(mlgm, sampleFreq=1)		

#plot(sData[[1]])

#RUs<-sData[[1]]@associationData[,2]
#Ts<-sData[[1]]@associationData[,1]
Ts<-seq(0,100,by=0.1)
#start testing the polynomial 
kon<-c(3e3)
koff<-c(5E-2)#,2E-2)
Rmax<-c(60)
conc<-1E-4
degree<-110
associationLength<-100 
dissociationLength<-100

#first testing transform ts with factorial
tf<-preparePolynomialTsWithFactorial(Ts, degree);

#get polynomial coefficient
spc<-getSprPolynomialCoef(kon,koff, Rmax, conc, degree)

#get SPR polynomial approximated
paSPR<-getPolynomialApproximatedSPR(kon, koff, Rmax,Ts, conc, degree)

#now start doing simulation with multiligand model, with a single component
mlgm<- new("MultiLigandModel", kon=kon, koff=koff, analyteConcentrations=conc, 
		associationLength=associationLength, dissociationLength=dissociationLength, Rmax=Rmax
		)
sData<-Simulate(mlgm, sampleFreq=1)	
Ts<-sData[[1]]@associationData[,1]
tf<-preparePolynomialTsWithFactorial(Ts, degree);
paSPR<-getPolynomialApproximatedSPR(kon, koff, Rmax,Ts, conc, degree)

#plot them together
 plot(sData[[1]])

lines(Ts,paSPR, col=2, lty=2)
legend(50,40, c("simulation by multiLigand model","polynomial with degree of 110"), col=c(1,2), pch=c(-1,-1), lty=c(1,2))


sData[[1]]@associationData
sfitPoly<-lm(formula = sData[[1]]@associationData[,2] ~ tf +0
	, weights=exp(-1*Ts)
	#, weights = 1/factorial(floor(x + 0.1))
	)
	#fitting
fpc<-fitPolySPR(Ts, sData[[1]]@associationData[,2])
predictFitSPR<-predictPolynomialApproximatedSPR(Ts,fpc[["lm.coefficient"]])
lines(Ts, predictFitSPR, col=3, lty=3)
plotFitSPR.kon(sData[[1]], fpc)

####--------->
#testing fitPolySPRs
#for non-steady state spr
kon<-c(2E3, 1E3, 3E3)
koff<-c(0.001, 0.01,0.05)

analyteConcentrations<-c(
							#5E-4, 4E-4, 3E-4, 2E-4, 
							1E-4, 9E-5,6E-5, 3E-5, 2E-5)
associationLength<-2000 
dissociationLength<-2000
Rmax<-c(50,80,60)

mlgm<- new("MultiLigandModel", kon=kon, koff=koff, analyteConcentrations=analyteConcentrations, 
		associationLength=associationLength, dissociationLength=dissociationLength, Rmax=Rmax
		)
sData<-Simulate(mlgm, sampleFreq=1)	
fpc<-fitPolySPRs(sData[[1]], debug=T, weights.step=1, weights.scale=1);

spc<-getSprPolynomialCoef(kon,koff, Rmax, analyteConcentrations[1], degree=100)

###======>
###following the above code, starting testing the fitting coefficient of polySPR
#first get the distribution of Rmax
E_kon<-sum(Rmax/sum(Rmax)*kon)
E_kon2<-sum(Rmax/sum(Rmax)*kon^2)
E_kon3<-sum(Rmax/sum(Rmax)*kon^3)
E_kon4<-sum(Rmax/sum(Rmax)*kon^4)
E_kon5<-sum(Rmax/sum(Rmax)*kon^5)

m<-fitCoefficientPolySPRs(fpc, analyteConcentrations)
m/sum(Rmax) ##this is the estimated, since m is Rmax scaled distribution

###=======> now following the previous code, we start testing the koffs
m1<-fitPolySPRsKoff(sData[[1]], debug=T, weights.type="factorial", degree=200, weights.step=10)

m2<-fitPolySPRsKoff(sData[[1]], debug=T, weights.type="exp", degree=200)

m3<-fitPolySPRsKoff(sData[[1]], debug=T, weights.type="uniform", degree=200)

#start to get the distribution from koffs
r0<-matrix(Rmax,nrow=length(analyteConcentrations),ncol=length(Rmax), byrow=T)
for(i in 1:length(analyteConcentrations))
{
r0[i,]<-Rmax*kon*analyteConcentrations[i]/(kon*analyteConcentrations[i]+koff)*(1-exp(-(kon*analyteConcentrations[i]+koff)*associationLength))
}
#expected koff moments
E_R0<-c(sum(r0[1,]),sum(r0[2,]),sum(r0[3,]),sum(r0[4,]),sum(r0[5,]))
E_koff1<-rep(0,length(analyteConcentrations))
E_koff2<-rep(0,length(analyteConcentrations))
E_koff3<-rep(0,length(analyteConcentrations))
E_koff4<-rep(0,length(analyteConcentrations))

for(i in 1:length(analyteConcentrations))
{
	E_koff1[i]<-sum(koff*r0[i,]/sum(r0[i,]))
	E_koff2[i]<-sum(koff^2*r0[i,]/sum(r0[i,]))
	E_koff3[i]<-sum(koff^3*r0[i,]/sum(r0[i,]))
	E_koff4[i]<-sum(koff^4*r0[i,]/sum(r0[i,]))
	
}
#fitted koff moments
f_E_koff1<-rep(0,length(analyteConcentrations))
f_E_koff2<-rep(0,length(analyteConcentrations))
f_E_koff3<-rep(0,length(analyteConcentrations))
f_E_koff4<-rep(0,length(analyteConcentrations))
for(i in 1:length(analyteConcentrations))
{
	f_E_koff1[i]<-m[[i]][[1]][2]#/m[[i]][[1]][1]
	f_E_koff2[i]<-m[[i]][[1]][3]#/m[[i]][[1]][1]
	f_E_koff3[i]<-m[[i]][[1]][4]#/m[[i]][[1]][1]
	f_E_koff4[i]<-m[[i]][[1]][5]#/m[[i]][[1]][1]
	
}

###now testing the steady state koff
##based on the previous fitPolySPRsKoff
mss<-fitCoefficientPolySPRsKoff(m, analyteConcentrations, sum(Rmax))

#expected koff moments
E_koff<-rep(0,7)
E_koff[1]<-sum(Rmax)
E_koff[2]<-sum(koff*Rmax/sum(Rmax))
E_koff[3]<-sum(koff^2*Rmax/sum(Rmax))
E_koff[4]<-sum(koff^3*Rmax/sum(Rmax))
E_koff[5]<-sum(koff^4*Rmax/sum(Rmax))
E_koff[6]<-sum(koff^5*Rmax/sum(Rmax))
E_koff[7]<-sum(koff^6*Rmax/sum(Rmax))