Core of code for Anatomic set of mechanical models for the organ of Corti

Posted in bioRxiv https://doi.org/10.1101/760835

For further details, you are welcomed to contact the authors.

Geometric parameters

L=10;(*lengthofRL*)LH=10;(*lengthofHC*)LT=L+0.5LH^2/(L+LH);(*lengthofTM*)x1=3;(*centerof1stOHC*)x2=7;(*centerof2ndOHC*)ell=2;(*extensionofOHCinradialdirection*)

Dynamic parameters

m=10;(*massateachendofanOHC*)mH=120;(*massofHC*)IRL=2.*^3;(*momentofinertiaofRL*)KRL=1.*^3;(*torsionalspringconstantonRL*)KCP=1.*^2;(*springconstantbetweenRLandCP*)betaCP=3;(*dampingconstantbetweenRLandCP*)KC=50;(*springconstantbetweenendsofOHC*)betaC=43;(*dampingconstantbetweenendsofOHC*)DEL=0.25375;(*maximumcontraction*)KD1=4.*^2;(*springconstantbetween1stDCandBM*)KD2=4.*^2;(*springconstantbetween2ndDCandBM*)beta1=3;(*damping1stDC*)beta2=3;(*damping2ndDC*)beta12=3;(*frictionbetweenDC1-2*)KB=10;(*springconstantouterbundles*)H1=6.5*^-3;(*extensionunstablerangeOHB1*)H2=5.*^-3;(*extensionunstablerangeOHB2*)KIHC=10;(*torsionalspringconstantinnerbundle*)HIN=5.*^-3;(*extensionsoftrangeIHB*)cs=2;(*innersulcusresistancetocompression*)mu=10;(*dragcoefficientIHB*)amp=310^-5;(*inputamplitude*)wbm=5.3294;(*inputangularfrequency*)noise=2.*^-6;(*representativevaluefor2k_BT/thickness*)

Variables

q: Q(0, t); tet: θ; b1:
b
1
; b2:
b
2
; s1:
s
1
; s2:
s
2
; tin:
θ
in
; pin:
p
in

wn1, wn2, wn3, wn4: noise frequencies
ph1, ph2, ph3, ph4: noise phases

Initial values

(*Letter"p"standsfortimederivativeand"0"forinitialvalue*)tet0=RandomVariate[NormalDistribution[0,Sqrt[noise/KRL]]];tetp0=RandomVariate[NormalDistribution[0,Sqrt[noise/IRL]]];b10=RandomVariate[NormalDistribution[0,Sqrt[noise/KCP]]];b1p0=RandomVariate[NormalDistribution[0,Sqrt[noise/m]]];b20=RandomVariate[NormalDistribution[0,Sqrt[noise/KCP]]];b2p0=RandomVariate[NormalDistribution[0,Sqrt[noise/m]]];s10=RandomVariate[NormalDistribution[0,Sqrt[noise/KD1]]];s1p0=RandomVariate[NormalDistribution[0,Sqrt[noise/m]]];s20=RandomVariate[NormalDistribution[0,Sqrt[noise/KD1]]];s2p0=RandomVariate[NormalDistribution[0,Sqrt[noise/m]]];pin0=RandomVariate[NormalDistribution[0,Sqrt[noisecs]]];q0=0;tin0=0;

Input

ybm=ampCos[wbmt];(*signal;caseofsteadyamplitude*)an=Sqrt[0.5noise/(KD1+KD2)];wn10=RandomReal[{0,2wbm}];wn20=RandomReal[{0,2wbm}];wn30=RandomReal[{0,2wbm}];wn40=RandomReal[{0,2wbm}];ph10=RandomReal[{0,2Pi}];ph20=RandomReal[{0,2Pi}];ph30=RandomReal[{0,2Pi}];ph40=RandomReal[{0,2Pi}];yn1=anCos[wn1[t]t-ph1[t]];yn2=anCos[wn2[t]t-ph2[t]];yn3=anCos[wn3[t]t-ph3[t]];yn4=anCos[wn4[t]t-ph4[t]];ytot=ybm+yn1+yn2+yn3+yn4;(*signal+noise*)

x-dependence, forces and torques

(*evaluationofthiscellmaytakeafewseconds*)ylow=tet[t]x+Piecewise[{{b1[t](1+Cos[2Pi(x-x1)/ell]),x1-ell/2<x<x1+ell/2},{b2[t](1+Cos[2Pi(x-x2)/ell]),x2-ell/2<x<x2+ell/2},{-tet[t](L+LH)(x-L)^2/LH^2,x>L}},0];qx=Integrate[Evaluate[D[ylow,t]],{x,0,xx},Assumptionsxx>0]/.xxx;(*Q(x,t)-Q(0,t)*)px=-12Integrate[Evaluate[qx+D[qx,t]/10]+q[t]+q'[t]/10,{x,LT,xx},Assumptionsxx<LT]/.xxx//Simplify;(*p(x,t)*)p0=px/.x0;torqueRL=-Integrate[pxx,{x,0,x1-ell/2}]-Integrate[pxx,{x,x1+ell/2,x2-ell/2}]-Integrate[pxx,{x,x2+ell/2,L}]//Simplify;(*torque,pressureonRL*)FH=-mH(5L+LH)tet''[t]/12-Integrate[pxx,{x,L,LT}]/LH//Simplify;(*force,HConRL*)Fp1=-Integrate[px,{x,x1-ell/2,x1+ell/2}];(*force,pressureonCP1*)Fp2=-Integrate[px,{x,x2-ell/2,x2+ell/2}];(*force,pressureonCP2*)F1=KCPb1[t]+betaCPb1'[t];(*force,CP1onRL*)F2=KCPb2[t]+betaCPb2'[t];(*force,CP2onRL*)bund=KBPiecewise[{{-h-H,h<-H},{-h+H,h>H}},(H/Pi)Sin[Pih/H]];(*auxiliary*)FB1=bund/.htet[t]x1+b1[t]/.HH1;(*force,bundleonCP1*)FB2=bund/.htet[t]x2+b2[t]/.HH2;(*force,bundleonCP1*)FC1=KC(tet[t]x1+b1[t]-s1[t]+DELTanh[(tet[t]x1+b1[t])/H1])+betaC(tet'[t]x1+b1'[t]-s1'[t]);(*force,OHConCP1*)FC2=KC(tet[t]x2+b2[t]-s2[t]+DELTanh[(tet[t]x2+b2[t])/H2])+betaC(tet'[t]x2+b2'[t]-s2'[t]);(*force,OHConCP2*)tauIHC=-KIHCPiecewise[{{tin[t]+HIN,tin[t]<-1.5HIN},{tin[t]-HIN,tin[t]>1.5HIN}},4tin[t]^3/(27HIN^2)];(*torque,IHConIHB*)

Dynamic equations

eqRL=Simplify[IRLtet''[t]-KRLtet[t]+F1x1+F2x2+FHL+torqueRL];eqCP1=Simplify[m(tet''[t]x1+b1''[t])-F1+FB1-FC1+Fp1];eqCP2=Simplify[m(tet''[t]x2+b2''[t])-F2+FB2-FC2+Fp2];eqD1=Simplify[ms1''[t]FC1+KD1(ytot-s1[t])-beta1s1'[t]+beta12(s2'[t]-s1'[t])];eqD2=Simplify[ms2''[t]FC2+KD2(ytot-s2[t])-beta2s2'[t]-beta12(s2'[t]-s1'[t])];eqIS=pin'[t]-csq[t];eqP=pin[t]-p0mu(q[t]+tin'[t]/2);eqIHB=tauIHCmu(q[t]/2+tin'[t]/3);

Time integration

tend=1000;(*desiredfollowuplapseoftime*)sol=NDSolve[{eqRL,eqCP1,eqCP2,eqD1,eqD2,eqIS,eqP,eqIHB,tet[0]tet0,tet'[0]tetp0,b1[0]b10,b1'[0]b1p0,b2[0]b20,b2'[0]b2p0,s1[0]s10,s1'[0]s1p0,s2[0]s20,s2'[0]s2p0,tin[0]tin0,q[0]q0,pin[0]pin0,wn1[0]wn10,wn2[0]wn20,wn3[0]wn30,wn4[0]wn40,ph1[0]ph10,ph2[0]ph20,ph3[0]ph30,ph4[0]ph40,WhenEvent[Mod[t,0.7]0.0,{prev=wn1[t],wn1[t]RandomReal[{0,2wbm}],ph1[t]Mod[ph1[t]+(wn1[t]-prev)t,2Pi]}],WhenEvent[Mod[t,0.9]0.07,{prev=wn2[t],wn2[t]RandomReal[{0,2wbm}],ph2[t]Mod[ph2[t]+(wn2[t]-prev)t,2Pi]}],WhenEvent[Mod[t,1.1]0.0707,{prev=wn3[t],wn3[t]RandomReal[{0,2wbm}],ph3[t]Mod[ph3[t]+(wn3[t]-prev)t,2Pi]}],WhenEvent[Mod[t,1.3]0.070707,{prev=wn4[t],wn4[t]RandomReal[{0,2wbm}],ph4[t]Mod[ph4[t]+(wn4[t]-prev)t,2Pi]}]},{tet,q,b1,b2,s1,s2,pin,tin,wn1,wn2,wn3,wn4,ph1,ph2,ph3,ph4},{t,0,tend},DiscreteVariables{wn1,wn2,wn3,wn4,ph1,ph2,ph3,ph4}][[1]];(*outputinformation*)

Re-initialization for extension of follow up time

tet0=tet[tend]/.sol;tetp0=tet'[tend]/.sol;b10=b1[tend]/.sol;b1p0=b1'[tend]/.sol;b20=b2[tend]/.sol;b2p0=b2'[tend]/.sol;s10=s1[tend]/.sol;s1p0=s1'[tend]/.sol;s20=s2[tend]/.sol;s2p0=s2'[tend]/.sol;tin0=tin[tend]/.sol;pin0=pin[tend]/.sol;q0=q[tend]/.sol;

Amplitude and phase

tinit=tend-7;(*beginningofdesiredlapseoftime*)tfin=tend-2.1Pi/wbm;(*endofdesiredlapseoftime*)func=tin[t]/.sol;(*pickdesiredfunction*)t1=t/.FindMaximum[func,{t,tinit,tinit+0.1/wbm}][[2]];(*refinelapseoftimetoobtainintegernumberofperiods*)t2=t/.FindMaximum[func,{t,tfin,tfin-0.1/wbm}][[2]];average=NIntegrate[func,{t,t1,t2},Method"Oscillatory",MaxRecursion20,PrecisionGoal2]/(t2-t1);averagesq=NIntegrate[func^2,{t,t1,t2},Method"Oscillatory",MaxRecursion20]/(t2-t1);amplitudesq=averagesq-average^2;(*thisisthesquareoftheamplitude*)phase=ArcTan[NIntegrate[funcCos[wbmt],{t,t1,t2}],-NIntegrate[funcSin[wbmt],{t,t1,t2}]];(*phaseoffunc*)

Cite this as: Jorge Berger, "Models for Organ of Corti" from the Notebook Archive (2020), https://notebookarchive.org/2020-01-2aqmevw

Wolfram Foundation

Download