From 11c2855f4b82bd3eae86535430722914a46e28ba Mon Sep 17 00:00:00 2001 From: "Cody, Jonathan William" Date: Sun, 25 Dec 2022 08:42:47 -0600 Subject: [PATCH] Delete N803_varied_reg_fA.m --- ICsolvers/N803_varied_reg_fA.m | 283 --------------------------------- 1 file changed, 283 deletions(-) delete mode 100644 ICsolvers/N803_varied_reg_fA.m diff --git a/ICsolvers/N803_varied_reg_fA.m b/ICsolvers/N803_varied_reg_fA.m deleted file mode 100644 index cc63047..0000000 --- a/ICsolvers/N803_varied_reg_fA.m +++ /dev/null @@ -1,283 +0,0 @@ -%% N803_varied_reg_fA.m - solves model for 2 cohorts with shared reg and fA -% -% /--------------------------------------------------------------\ -% | Date: 07/01/2022 | -% | Author: Jonathan Cody | -% | Affiliation: Purdue University | -% | Weldon School of Biomedical Engineering | -% | Pienaar Computational Systems Pharmacology Lab | -% \--------------------------------------------------------------/ -% -% Nomenclature: V = SIV virions [#/無] -% T8 = total CD8+ T cells [#/無] -% E0 = resting SIV-specific CD8+ T cells [#/無] -% Ea = active SIV-specific CD8+ T cells [#/無] -% B0 = resting bystander CD8+ T cells [#/無] -% Ba = active bystander CD8+ T cells [#/無] -% X = N803 at absorption site [pmol/kg] -% C = N803 plasma concentration [pM] -% R = regulation [] (dimensionless quantity) -% -%% ======================================================================== -% INPUTS -% ======================================================================== -% -% SoluTimes = ascending vector of days at which to evaluate solution -% -% DoseTimes{c} = ascending vector of days at which to administer doses -% (elements of 'DoseTimes' must also be in 'SoluTimes') -% -% AllPars = vector of parameters (see list in function) -% -%% ======================================================================== -% OPTIONS -% ======================================================================== -% -% SkipTimes{c} = [min max] time point beyond which to skip model solving -% (outputs before 'min' will be made equal to output at 'min') -% (leave as [] to ignore and solve for all 'SoluTimes') -% -% oneCohort = scalar to run model for just one cohort ('1' or '2') -% (leave as [] to ignore and solve for both cohorts) -% -% All additional inputs will be passed as a cell vector to 'N803_model_2' -% and used to define options (see function for list) -% EX: N803_single(SoluTimes,DoseTimes,AllPars,'AbsTol',1e-2} -% will set ode solver absolute tolerance to 1e-2 -% -%% ======================================================================== -% OUTPUTS -% ======================================================================== -% -% Y_OUT(:,1) = V at points in 'SoluTimes' [log fold change] cohort 1 -% Y_OUT(:,2) = T8 at points in 'SoluTimes' [fold change] cohort 1 -% Y_OUT(:,3) = V at points in 'SoluTimes' [log fold change] cohort 2 -% Y_OUT(:,4) = T8 at points in 'SoluTimes' [fold change] cohort 2 -% -% PARS(1,:) = [ parameters , initial conditions ] for cohort 1 (see code) -% PARS(2,:) = [ parameters , initial conditions ] for cohort 2 (see code) -% -%% ======================================================================== -% FUNCTION -% ======================================================================== -function [Y_OUT,PARS] = ... - N803_varied_reg_fA(SoluTimes,DoseTimes,AllPars,SkipTimes,oneCohort,... - varargin) - -if isempty(oneCohort) ; RunCohort = [ 1 1 ] ; -elseif oneCohort == 1 ; RunCohort = [ 1 0 ] ; -else ; RunCohort = [ 0 1 ] ; -end - -Y_Cohort = cell(1,2) ;% cell for storing outputs -P_Cohort = cell(1,2) ;% cell for storing parmeters - -% Rename inputed parameters ----------------------------------------------- -Vi(1) = AllPars(01) ;% V initial value [log(#/mL)] (cohort 1) -Vi(2) = AllPars(02) ;% V initial value [log(#/mL)] (cohort 2) -EBi(1) = AllPars(03) ;% E+B initial value [#/無] (cohort 1) -EBi(2) = AllPars(04) ;% E+B initial value [#/無] (cohort 2) -fE(1) = AllPars(05) ;% initial frequency: E/(E+B) (cohort 1) -fE(2) = AllPars(06) ;% initial frequency: E/(E+B) (cohort 2) - -q = AllPars(07) ;% V growth rate (if E+B were absent) [/d] -psi = AllPars(08) ;% ratio of base killing rates gB0/gE0 -V50E = AllPars(09) ;% 50% viral stimulation saturation for E [#/mL] -V50B = AllPars(10) ;% 50% viral stimulation saturation for B [#/mL] -mEn = AllPars(11) ;% normalized Ea reversion rate constant [] -mBn = AllPars(12) ;% normalized Ba reversion rate constant [] - -EB50 = AllPars(13) ;% 50% E+B proliferation saturation [#/無] -p0n = AllPars(14) ;% nomalized E0/B0 proliferation rate constant [/d] -pE = AllPars(15) ;% Ea proliferation rate constant [/d] -pB = AllPars(16) ;% Ba proliferation rate constant [/d] -d = AllPars(17) ;% E0/B0 death rate constant [/d] -dA = AllPars(18) ;% Ea/Ba death rate constant [/d] - -Xi = AllPars(19) ;% X initial condition [pmol/kg] -ka = AllPars(20) ;% N803 absorption rate constant [/d] -ke = AllPars(21) ;% N803 elimination rate constant [/d] -vd = AllPars(22) ;% N803 'volume of distribution'/'bioavailability' [L/kg] -C50 = AllPars(23) ;% 50% N803 stimulation concentration [pM] (Cohort 1) -pm = AllPars(24) ;% E0/B0 maximum proliferation rate [] -aE1 = AllPars(25) ;% E activation stimulation factor [] -aB1 = AllPars(26) ;% B activation stimulation factor [] - -dR(1) = AllPars(27) ;% R decay rate constant [/d] (cohort 1) -dR(2) = AllPars(28) ;% R decay rate constant [/d] (cohort 2) -sig(1) = AllPars(29) ;% ratio of initial regulation due to B/E (cohort 1) -sig(2) = AllPars(30) ;% ratio of initial regulation due to B/E (cohort 2) -p2(1) = AllPars(31) ;% E0/B0 proliferation regulation factor [] (cohort 1) -p2(2) = AllPars(32) ;% E0/B0 proliferation regulation factor [] (cohort 2) -aE2s = AllPars(33) ;% E activation regulation factor [] (sampled) -aB2s = AllPars(34) ;% B activation regulation factor [] (sampled) -gE2s = AllPars(35) ;% E killing regulation factor [] (sampled) -gB2(1) = AllPars(36) ;% B killing regulation factor [] (cohort 1) -gB2(2) = AllPars(37) ;% B killing regulation factor [] (cohort 2) - -%% ------------------------------------------------------------------------ -% Calculate some initial conditions & parameters -------------------------- - -Vi = 10.^(Vi - 3) ;% converting V initial value [#/無] -V50E = V50E/1000 ;% 50% viral stimulation saturation for E [#/無] -V50B = V50B/1000 ;% 50% viral stimulation saturation for B [#/無] -YE = Vi ./ (Vi+ V50E) ;% initial V/(V50E+V) [cohort 1,2] -YB = Vi ./ (Vi+ V50B) ;% initial V/(V50B+V) [cohort 1,2] - -% restrict mE and mB such that initial activation aE and aB are positive -UE = (2*pE/(pE+dA))^7 ; -UB = 2*pB/(pB+dA) ; -mE = mEn*dA/(UE-1) ;% Ea reversion rate constant [/d] -mB = mBn*dA/(UB-1) ;% Ba reversion rate constant [/d] - -% solve for initial ratios below (based on active steady-state) -ZE = UE/(mE+dA) ; -for i = 1:7 - ZE = ZE + (2*pE)^(i-1)/(pE+dA)^i ;% EAi/aEi/E0i -end -ZB = 1/(pB+dA) + UB/(mB+dA) ;% BAi/aBi/B0i - -WE = 1 ; -for i = 1:7 - WE = WE + (mE+dA)*(pE+dA)^(i-1)/(2*pE)^i ;% EAi/E8i -end -WB = 1 + (mB+dA)/(2*pB) ;% BAi/B2i - -QE = mE/WE - 1/ZE ;% collection -QB = mB/WB - 1/ZB ;% collection - -% restrict p0 to ensure same sign for E0/EA and for B0/BA -p0_max = min( d*(1+p2).*(EB50+EBi)/EB50 ) ;% maximum possible value for p0 -p0 = p0n * p0_max ;% E0/B0 prolif rate constant [/d] -p1 = pm/p0 ;% E0/B0 proliferation stimulation factor -pi = p0*EB50./(EB50+EBi)./(1+p2) ;% initial E0/B0 proliferation rate [/d] - -% calculate aE0,aB0 (same for each cohort) and aE2,aB2 (different) -aEi = (d-pi)/QE/ZE ;% initial E0 activation rate [/d] -aBi = (d-pi)/QB/ZB ;% initial E0 activation rate [/d] - -aE_YE_ratio = aEi ./ YE ; -if aE_YE_ratio(1) > aE_YE_ratio(2) ; s = 1 ; c = 2 ; -else ; s = 2 ; c = 1 ; -end -aE2(s) = aE2s ;% regulation of E activation (sampled) -aE2(c) = aE_YE_ratio(s) / aE_YE_ratio(c) * (1+aE2s) - 1 ;% (calculated) - -aB_YB_ratio = aBi ./ YB ; -if aB_YB_ratio(1) > aB_YB_ratio(2) ; s = 1 ; c = 2 ; -else ; s = 2 ; c = 1 ; -end -aB2(s) = aB2s ;% regulation of B activation (sampled) -aB2(c) = aB_YB_ratio(s) / aB_YB_ratio(c) * (1+aB2s) - 1 ;% (calculated) - -aE0 = aE_YE_ratio .* (1+aE2) ;% base activation constant for E -aB0 = aB_YB_ratio .* (1+aB2) ;% base activation constant for B - -% solve for initial values of E0,Ea,B0,Ba each cohort -Ei = EBi.*fE ;% initial E -Bi = EBi.*(1-fE) ;% initial B -E0 = Ei./(1+ZE*aEi) ;% initial E0 -B0 = Bi./(1+ZB*aBi) ;% initial B0 -EA = E0.*(ZE*aEi) ;% initial EA -BA = B0.*(ZB*aBi) ;% initial BA -Ri = 1 ;% regulation normalized to initial value - -% calculate gE0,gB0 (same for each cohort) and gE2 (different) -if EA(1) < EA(2) ; s = 1 ; c = 2 ; -else ; s = 2 ; c = 1 ; -end -gE2(s) = gE2s ;% regulation of E killing (sampled) -beta = psi .* BA ./ (1+gB2) ; -eps(s) = EA(s) / (1+gE2s) ; -eps(c) = eps(s) + beta(s) - beta(c) ; -gE2(c) = EA(c) / eps(c) - 1 ;% regulation of E killing (calculated) - -gE0 = q ./ ( eps + beta ) ;% base killing constant for E -gB0 = psi * gE0 ;% base killing constant for B - -% calculate sE,sB (different for each cohort) -sE = dR ./ ( YE + sig.*YB ) ;% reg gen by E act -sB = sig.*sE ;% regulation generation by B activation - -%% Do for each cohort (NOT indenting loop) ================================ -for c = 1:2 - -% solve for initial E1-8 and B1-2 -E = zeros(1,8) ;% initial E1-E8 -E(8) = EA(c)/WE ;% E8 -E(7) = E(8) * (mE+dA)/(2*pE) ;% E7 -for i = 6:-1:1 - E(i) = E(i+1) * (pE+dA) / (2*pE) ;% E6 to E1 -end -B = BA(c)/WB * [ (mB+dA) / (2*pB) , 1 ] ;% initial B1-B2 - -%% ------------------------------------------------------------------------ -% Prepare parameter and initial value vectors and call 'N803_model_2' ----- - -Pars(01) = q ;% V growth rate (if E+B were absent) [/d] -Pars(02) = gE0(1) ;% Ea killing rate constant [無/#-d] -Pars(03) = gB0(1) ;% Ba killing rate constant [無/#-d] - -Pars(04) = V50E ;% 50% viral stimulation saturation for E [#/無] -Pars(05) = V50B ;% 50% viral stimulation saturation for B [#/無] -Pars(06) = aE0(1) ;% E0 activation rate constant [/d] -Pars(07) = aB0(1) ;% B0 activation rate constant [/d] -Pars(08) = mE ;% Ea reversion rate constant [/d] -Pars(09) = mB ;% Ba reversion rate constant [/d] - -Pars(10) = EB50 ;% 50% E+B proliferation saturation [#/無] -Pars(11) = p0 ;% E0/B0 proliferation rate constant [/d] -Pars(12) = pE ;% Ea proliferation rate constant [/d] -Pars(13) = pB ;% Ba proliferation rate constant [/d] -Pars(14) = d ;% E0/B0 death rate constant [/d] -Pars(15) = dA ;% Ea/Ba death rate constant [/d] - -Pars(16) = Xi ;% X initial condition [pmol/kg] -Pars(17) = ka ;% N803 absorption rate constant [/d] -Pars(18) = ke ;% N803 elimination rate constant [/d] -Pars(19) = vd ;% N803 'volume of distribution'/'bioavailability' [L/kg] -Pars(20) = C50 ;% 50% N803 stimulation concentration [pM] -Pars(21) = p1 ;% E0/B0 proliferation stimulation factor [] -Pars(22) = aE1 ;% E activation stimulation factor [] -Pars(23) = aB1 ;% B activation stimulation factor [] - -Pars(24) = sE(c) ;% R generation due to E0 activation [/d] -Pars(25) = sB(c) ;% R generation due to B0 activation [/d] -Pars(26) = dR(c) ;% R decay rate constant [/d] -Pars(27) = gE2(c) ;% E killing regulation factor [] -Pars(28) = gB2(c) ;% B killing regulation factor [] -Pars(29) = p2(c) ;% E0/B0 proliferation regulation factor [] -Pars(30) = aE2(c) ;% E activation regulation factor [] -Pars(31) = aB2(c) ;% B activation regulation factor [] - -% V E0-8 B0-2 X C R initial values -Yic = [ Vi(c) E0(c) E B0(c) B 0 0 Ri Ri ] ; - -if any( [ Pars Yic ] < 0 ) - error('Negative parameters or initial values.') -end - -% If 'SkipTimes' is empty, do not skip any times -if isempty(SkipTimes{c}) ; SkipTimes{c} = [-inf inf] ; end - -idLo = SoluTimes < SkipTimes{c}(1) ;% index of early times to skip soln -idHi = SoluTimes > SkipTimes{c}(2) ;% index of later times to skip soln -idSol = ~ ( idLo | idHi ) ;% index of times in 'SoluTimes' to solve - -if RunCohort(c) == 1 - Y_COH = N803_model_2(SoluTimes(idSol),DoseTimes{c},Pars,Yic,varargin) ; - Y_LO = ones(sum(idLo),1)*Y_COH(1 ,:) ;% constant Y for early times - Y_HI = ones(sum(idHi),1)*Y_COH(end,:) ;% constant Y for later times - Y_COH = [ Y_LO ; Y_COH ; Y_HI ] ;%#ok % total 'solution' matrix - Y_Cohort{c} = Y_COH ; - P_Cohort{c} = [ Pars Yic ] ; -end - -end - -%% Prepare main outputs 'Y_OUT' and 'PARS' ================================ - -Y_OUT = [ Y_Cohort{1} , Y_Cohort{2} ] ; -PARS = [ P_Cohort{1} ; P_Cohort{2} ] ; - -end \ No newline at end of file