import numpy import functools def fval(t,W): if callable(W): return W(t) return W def L(t,t0,W): ft0=fval(t,t0) fW=fval(t,W) x=(t-ft0)/fW if x<-18: return 0 if x>18: return 1 return 1/(1+numpy.exp(-x)) def dL(t,t0,W): return L(t,t0,W)*(1-L(t,t0,W)) def dLdW(t,t0,W): fW=fval(t,W) ft0=fval(t,t0) return dL(t,t0,W)*(t-ft0)/fW/fW def dLdt0(t,t0,W): fW=fval(t,W) return dL(t,t0,W)/fW def logisticPar(par): t0=par['t0']['value'] W=par['W']['value'] f=functools.partial(L,t0=t0,W=W) #D1[x] could be a function !!! Dt0=par['t0']['derivatives'] dfdt0=functools.partial(dLdt0,t0=t0,W=W) dt0={x:lambda t,df=dfdt0,D=Dt0,x=x: df(t)*to(D[x])(t) for x in Dt0} #D2[x] could be a function !!! DW=par['W']['derivatives'] dfdW=functools.partial(dLdW,t0=t0,W=W) dW={x:lambda t,df=dfdW,D=DW,x=x: df(t)*to(D[x])(t) for x in DW} #merge return Object(f,[dt0,dW]) def deltaW(t,t1,t2): ft1=fval(t,t1) ft2=fval(t,t2) return ft2-ft1 def deltaT(t,t1): ft1=fval(t,t1) return t-t1 def fq(t,t1,t2): return deltaT(t,t1)/deltaW(t,t1,t2) def f1(t,C,t1,t2,L1,L2): fC=fval(t,C) return fC*fq(t,t1,t2)*L1(t)*(1-L2(t)) def f2(t,C,L1,L2): fC=fval(t,C) return fC*L1(t)*L2(t) def lg(t,C,t1,t2,eps,L1,L2): fC=fval(t,C) return eps*fC+f1(t,C,t1,t2,L1,L2)+f2(t,C,L1,L2) def dlgdt1(t,C,t1,t2,L1,L2): fC=fval(t,C) return -fC*L1(t)*(1-L2(t))*(1+fq(t,t1,t2))/deltaW(t,t1,t2) def dlgdt2(t,C,t1,t2,L1,L2): return -f1(t,C,t1,t2,L1,L2)/deltaW(t,t1,t2) def dlgdL1(t,C,t1,t2,L1,L2): fC=fval(t,C) fL2=L2(t) return fC*(fq(t,t1,t2)*(1-fL2)+fL2) def dlgdL2(t,C,t1,t2,L1,L2): fC=fval(t,C) return fC*L1(t)*(1-fq(t,t1,t2)) def dlgdC(t,C,t1,t2,L1,L2): fL2=L2(t) return L1(t)*(fq(t,t1,t2)*(1-fL2)+fL2) def lgFixedSlope(t,alpha,t1,t2,eps,L1,L2): fA=fval(t,alpha) ft1=fval(t,t1) dW=deltaW(t,t1,t2) return fA*(dW*eps+L1(t)*((t-ft1)*(1-L2(t))+dW*L2(t))) def dlgFixedSlopedt1(t,alpha,L1,L2): fA=fval(t,alpha) return -fA*L1(t) def dlgFixedSlopedt2(t,alpha,L1,L2): fA=fval(t,alpha) return -fA*L1(t)*L2(t) def dlgFixedSlopedL1(t,alpha,t1,t2,L1,L2): return lgFixedSlope(t,alpha,t1,t2,L1,L2)/L1(t) def dlgFixedSlopedL2(t,alpha,t1,t2,L1,L2): fA=fval(t,alpha) ft2=fval(t,t2) return fA*L1(t)*(ft2-t) def dlgFixedSlopedalpha(t,alpha,t1,t2,L1,L2): fA=fval(t,alpha) return lgFixedSlope(t,alpha,t1,t2,L1,L2)/fA def linearGrowth(par): #make sure we never return zero or negative C=par['C']['value'] t1=par['t1']['value'] t2=par['t2']['value'] eps=1e-8 pL1=logisticPar({'t0':par['t1'],'W':par['W1']}) pL2=logisticPar({'t0':par['t2'],'W':par['W2']}) L1=pL1['value'] L2=pL2['value'] f=functools.partial(lg,C=C,t1=t1,t2=t2,eps=eps,L1=L1,L2=L2) Dt1=par['t1']['derivatives'] dfdt1=functools.partial(dlgdt1,C=C,t1=t1,t2=t2,L1=L1,L2=L2) dt1={x:lambda t,df=dfdt1,D=Dt1,x=x:df(t)*to(D[x])(t) for x in Dt1} Dt2=par['t2']['derivatives'] dfdt2=functools.partial(dlgdt2,C=C,t1=t1,t2=t2,L1=L1,L2=L2) dt2={x:lambda t,df=dfdt2,D=Dt2,x=x:df(t)*to(D[x])(t) for x in Dt2} DL1=pL1['derivatives'] dfdL1=functools.partial(dlgdL1,C=C,t1=t1,t2=t2,L1=L1,L2=L2) dL1={x:lambda t,df=dfdL1,D=DL1,x=x:df(t)*D[x](t) for x in DL1} DL2=pL2['derivatives'] dfdL2=functools.partial(dlgdL2,C=C,t1=t1,t2=t2,L1=L1,L2=L2) dL2={x:lambda t,df=dfdL2,D=DL2,x=x: df(t)*D[x](t) for x in DL2} DC=par['C']['derivatives'] dfdC=functools.partial(dlgdC,C=C,t1=t1,t2=t2,L1=L1,L2=L2) dC={x:lambda t,df=dfdC,D=DC,x=x:df(t)*to(D[x])(t) for x in DC} return Object(f,[dt1,dt2,dL1,dL2,dC]) def linearGrowthFixedSlope(par): #make sure we never return zero or negative alpha=par['alpha']['value'] t1=par['t1']['value'] t2=par['t2']['value'] eps=1e-8 pL1=logisticPar({'t0':par['t1'],'W':par['W1']}) pL2=logisticPar({'t0':par['t2'],'W':par['W2']}) L1=pL1['value'] L2=pL2['value'] f=functools.partial(lgFixedSlope,alpha=alpha,t1=t1,t2=t2,eps=eps,L1=L1,L2=L2) Dt1=par['t1']['derivatives'] dfdt1=functools.partial(dlgFixedSlopedt1,alpha=alpha,L1=L1,L2=L2) dt1={x:lambda t,df=dfdt1,D=Dt1,x=x:df(t)*to(D[x])(t) for x in Dt1} Dt2=par['t2']['derivatives'] dfdt2=functools.partial(dlgFixedSlopedt2,alpha=alpha,L1=L1,L2=L2) dt2={x:lambda t,df=dfdt2,D=Dt2,x=x:df(t)*to(D[x])(t) for x in Dt2} DL1=pL1['derivatives'] dfdL1=functools.partial(dlgFixedSlopedL1,alpha=alpha,t1=t1,t2=t2,L1=L1,L2=L2) dL1={x:lambda t,df=dfdL1,D=DL1,x=x:df(t)*D[x](t) for x in DL1} DL2=pL2['derivatives'] dfdL2=functools.partial(dlgFixedSlopedL2,alpha=alpha,t1=t1,t2=t2,L1=L1,L2=L2) dL2={x:lambda t,df=dfdL2,D=DL2,x=x: df(t)*D[x](t) for x in DL2} DA=par['alpha']['derivatives'] dfdA=functools.partial(dlgFixedSlopedalpha,alpha=alpha,t1=t1,t2=t2,L1=L1,L2=L2) dA={x:lambda t,df=dfdA,D=DA,x=x:df(t)*to(D[x])(t) for x in DA} return Object(f,[dt1,dt2,dL1,dL2,dA]) def Object(v,dArray): o={'value':v,'function':True} allKeys=[] for x in dArray: allKeys.extend(x.keys()) allKeys=list(set(allKeys)) o['derivatives']={x:sumDict(dArray,x) for x in allKeys} return o def sumDict(dArray,x): #print('sumFunctions for {}'.format(x)) fs=[] for a in dArray: #print('\tChecking {}'.format(a)) try: #check if a[x] is a function with a dummy parameter fs.append(a[x]) #print('\tAdding [{}]'.format(v)) except KeyError: #print('\t{} not in table'.format(x)) continue return sumArray(fs) def sumArray(fs): return lambda t,fs=fs: sum([to(f)(t) for f in fs]) def to(a): if callable(a): return a return lambda t,a=a:a def contains(arr): for a in arr: if callable(a): return True return False def multiply(t,a,b): #abstract actual type of a and b fa=fval(t,a) fb=fval(t,b) return fa*fb def generate(df,D): return {x:functools.partial(multiply,a=D[x],b=df) for x in D}