From 7033560e2768129d5e2709e787ecf141170f4cf9 Mon Sep 17 00:00:00 2001 From: "Cody, Jonathan William" Date: Sun, 8 Jan 2023 14:42:54 -0600 Subject: [PATCH] Add files via upload --- N803_model_2.m | 109 ++++++++++++------------ N803_plotter_2_parameters.m | 114 ++++++++++++++++++++++++++ N803_shared.m | 42 +++++----- N803_single.m | 40 ++++----- main_shared_calibration_2022_08_18.m | 6 +- main_singles_calibration_2022_09_25.m | 4 +- 6 files changed, 214 insertions(+), 101 deletions(-) create mode 100644 N803_plotter_2_parameters.m diff --git a/N803_model_2.m b/N803_model_2.m index 28d06b2..63a29fc 100644 --- a/N803_model_2.m +++ b/N803_model_2.m @@ -1,19 +1,19 @@ %% N803_model_2.m - solves ODE model for N-803 treatment of SIV % % /--------------------------------------------------------------\ -% | Date: 12/25/2022 | +% | Date: 01/08/2022 | % | Author: Jonathan Cody | % | Affiliation: Purdue University | % | Weldon School of Biomedical Engineering | % | Pienaar Computational Systems Pharmacology Lab | % \--------------------------------------------------------------/ % -% Nomenclature: V = SIV virions [#/µL] -% T8 = total CD8+ T cells [#/µL] -% S0 = resting SIV-specific CD8+ T cells [#/µL] -% Sa = active SIV-specific CD8+ T cells [#/µL] -% N0 = resting non-SIV-specific CD8+ T cells [#/µL] -% Na = active non-SIV-specific CD8+ T cells [#/µL] +% Nomenclature: V = SIV virions [#/uL] +% T8 = total CD8+ T cells [#/uL] +% S0 = resting SIV-specific CD8+ T cells [#/uL] +% Sa = active SIV-specific CD8+ T cells [#/uL] +% N0 = resting non-SIV-specific CD8+ T cells [#/uL] +% Na = active non-SIV-specific CD8+ T cells [#/uL] % X = N803 at absorption site [pmol/kg] % C = N803 plasma concentration [pM] % R = regulation [] (dimensionless quantity) @@ -94,17 +94,17 @@ % Rename inputed parameters ----------------------------------------------- q = Pars(01) ;% V growth rate (if S+N were absent) [/d] -gS0 = Pars(02) ;% Sa killing rate constant [µL/#-d] -gN0 = Pars(03) ;% Na killing rate constant [µL/#-d] +gS0 = Pars(02) ;% Sa killing rate constant [uL/#-d] +gN0 = Pars(03) ;% Na killing rate constant [uL/#-d] -V50S = Pars(04) ;% 50% viral stimulation saturation for S [#/µL] -V50N = Pars(05) ;% 50% viral stimulation saturation for N [#/µL] +V50S = Pars(04) ;% 50% viral stimulation saturation for S [#/uL] +V50N = Pars(05) ;% 50% viral stimulation saturation for N [#/uL] aS0 = Pars(06) ;% S0 activation rate constant [/d] aN0 = Pars(07) ;% N0 activation rate constant [/d] mS = Pars(08) ;% Sa reversion rate constant [/d] mN = Pars(09) ;% Na reversion rate constant [/d] -SN50 = Pars(10) ;% 50% S+N proliferation saturation [#/µL] +SN50 = Pars(10) ;% 50% S+N proliferation saturation [#/uL] p0 = Pars(11) ;% S0/N0 proliferation rate constant [/d] pS = Pars(12) ;% Sa proliferation rate constant [/d] pN = Pars(13) ;% Na proliferation rate constant [/d] @@ -184,11 +184,11 @@ Y_OUT = [ Y_OUT, Y, zeros(Pieces(end),46) ] ; -V = Y(:,01) ;% SIV virions [#/µL] -S0 = Y(:,02) ;% resting specific CD8+ T cells [#/µL] -N0 = Y(:,11) ;% resting non-SIV-specific CD8+ T cells [#/µL] -Sa = sum( Y(:,03:10) ,2) ;% active specific CD8+ T cells [#/µL] -Na = sum( Y(:,12:13) ,2) ;% active non-SIV-specific T cells [#/µL] +V = Y(:,01) ;% SIV virions [#/uL] +S0 = Y(:,02) ;% resting specific CD8+ T cells [#/uL] +N0 = Y(:,11) ;% resting non-SIV-specific CD8+ T cells [#/uL] +Sa = sum( Y(:,03:10) ,2) ;% active specific CD8+ T cells [#/uL] +Na = sum( Y(:,12:13) ,2) ;% active non-SIV-specific T cells [#/uL] Saf = Y(:,10) ;% last generation of Sa Naf = Y(:,13) ;% last generation of Na Sap = Sa - Saf ;% Sa that is proliferating @@ -198,12 +198,11 @@ p = F(:,13) ;% S0/N0 proliferation rate [/d] aS = F(:,14) ;% S0 activation rate [/d] aN = F(:,15) ;% N0 activation rate [/d] -gS = F(:,16) ;% Sa killing rate [#/µL-d] -gN = F(:,17) ;% Na killing rate [#/µL-d] +gS = F(:,16) ;% Sa killing rate [#/uL-d] +gN = F(:,17) ;% Na killing rate [#/uL-d] Y_OUT(:,20) = F(:,01) ;% N803 effect [] - -Y_OUT(:,21) = (S0+N0+Sa+Na) ;% CD8+ T cells [#/µL] +Y_OUT(:,21) = (S0+N0+Sa+Na) ;% CD8+ T cells [#/uL] Y_OUT(:,22) = Sa ;% Sa Y_OUT(:,23) = Na ;% Na Y_OUT(:,24) = (S0+Sa) ;% S @@ -214,35 +213,35 @@ Y_OUT(:,29) = log10( Sa ) ;% Sa (log) Y_OUT(:,30) = log10( Na ) ;% Na (log) -Y_OUT(:,31) = p.*S0 ;% S0 proliferation term [#/µL/d] -Y_OUT(:,32) = p.*N0 ;% N0 proliferation term [#/µL/d] -Y_OUT(:,33) = d.*S0 ;% S0 death term [#/µL/d] -Y_OUT(:,34) = d.*N0 ;% N0 death term [#/µL/d] +Y_OUT(:,31) = p.*S0 ;% S0 proliferation term [#/uL/d] +Y_OUT(:,32) = p.*N0 ;% N0 proliferation term [#/uL/d] +Y_OUT(:,33) = d.*S0 ;% S0 death term [#/uL/d] +Y_OUT(:,34) = d.*N0 ;% N0 death term [#/uL/d] Y_OUT(:,35) = p ;% S0/N0 proliferation rate [/d] Y_OUT(:,36) = F(:,02) ;% C multiplier for p* (log) Y_OUT(:,37) = F(:,05) ;% R multiplier for p* (log) Y_OUT(:,38) = F(:,10) ;% density-dep multiplier for p* (log) -Y_OUT(:,39) = aS.*S0 ;% S0 activation term [#/µL/d] -Y_OUT(:,40) = aN.*N0 ;% N0 activation term [#/µL/d] -Y_OUT(:,41) = pS*Sap ;% Sa proliferation term [#/µL/d] -Y_OUT(:,42) = pN*Nap ;% Na proliferation term [#/µL/d] -Y_OUT(:,43) = dA*Sa ;% Sa death term [#/µL/d] -Y_OUT(:,44) = dA*Na ;% Na death term [#/µL/d] -Y_OUT(:,45) = mS*Saf ;% Sa reversion term [#/µL/d] -Y_OUT(:,46) = mN*Naf ;% Na reversion term [#/µL/d] +Y_OUT(:,39) = aS.*S0 ;% S0 activation term [#/uL/d] +Y_OUT(:,40) = aN.*N0 ;% N0 activation term [#/uL/d] +Y_OUT(:,41) = pS*Sap ;% Sa proliferation term [#/uL/d] +Y_OUT(:,42) = pN*Nap ;% Na proliferation term [#/uL/d] +Y_OUT(:,43) = dA*Sa ;% Sa death term [#/uL/d] +Y_OUT(:,44) = dA*Na ;% Na death term [#/uL/d] +Y_OUT(:,45) = mS*Saf ;% Sa reversion term [#/uL/d] +Y_OUT(:,46) = mN*Naf ;% Na reversion term [#/uL/d] Y_OUT(:,47) = aS ;% S0 activation rate [/d] Y_OUT(:,48) = aN ;% N0 activation rate [/d] Y_OUT(:,49) = F(:,03) ;% C multiplier for aS* -Y_OUT(:,50) = ( F(:,04) + F(:,12) )./F(:,12) ;% C multiplier for aS* +Y_OUT(:,50) = ( F(:,04) + F(:,12) )./F(:,12) ;% C multiplier for aN* Y_OUT(:,51) = F(:,06) ;% R multiplier for aS* Y_OUT(:,52) = F(:,07) ;% R multiplier for aN* Y_OUT(:,53) = F(:,11) ;% V multiplier for aS* Y_OUT(:,54) = F(:,12) ;% V multiplier for aN* -Y_OUT(:,55) = q*V ;% V growth term [#/µL/d] -Y_OUT(:,56) = gS.*Sa.*V ;% S killing term [#/µL/d] -Y_OUT(:,57) = gN.*Na.*V ;% N killing term [#/µL/d] +Y_OUT(:,55) = q*V ;% V growth term [#/uL/d] +Y_OUT(:,56) = gS.*Sa.*V ;% S killing term [#/uL/d] +Y_OUT(:,57) = gN.*Na.*V ;% N killing term [#/uL/d] Y_OUT(:,58) = gS.*Sa ;% S killing rate [/d] Y_OUT(:,59) = gN.*Na ;% N killing rate [/d] Y_OUT(:,60) = F(:,08) ;% R multiplier for gS* @@ -270,11 +269,11 @@ error('ODE solver stalled.') end - V = Y(01) ;% SIV virions [#/µL] - S0 = Y(02) ;% resting specific CD8+ T cells [#/µL] - Sa = Y(03:10) ;% active specific CD8+ T cells [#/µL] (S1 to S8) - N0 = Y(11) ;% resting non-SIV-specific CD8+ T cells [#/µL] - Na = Y(12:13) ;% active non-SIV-specific T cells [#/µL] (N1 to N2) + V = Y(01) ;% SIV virions [#/uL] + S0 = Y(02) ;% resting specific CD8+ T cells [#/uL] + Sa = Y(03:10) ;% active specific CD8+ T cells [#/uL] (S1 to S8) + N0 = Y(11) ;% resting non-SIV-specific CD8+ T cells [#/uL] + Na = Y(12:13) ;% active non-SIV-specific T cells [#/uL] (N1 to N2) X = Y(14) ;% N803 at absorption site [pmol/kg] C = Y(15) ;% N803 plasma concentration [pM] R = Y(16:17) ;% immune regulation [] @@ -283,8 +282,8 @@ p = F(13) ;% S0/N0 proliferation rate [/d] aS = F(14) ;% S0 activation rate [/d] aN = F(15) ;% N0 activation rate [/d] - gS = F(16) ;% Sa killing rate [#/µL-d] - gN = F(17) ;% Na killing rate [#/µL-d] + gS = F(16) ;% Sa killing rate [#/uL-d] + gN = F(17) ;% Na killing rate [#/uL-d] s = F(18)+F(19) ;% regulation generation rate [/d] dY = Y ;% dY/dt @@ -307,19 +306,19 @@ function J = jacob(~,Y) - V = Y(01) ;% SIV virions [#/µL] - S0 = Y(02) ;% resting specific CD8+ T cells [#/µL] - Sa = sum(Y(03:10)) ;% active specific CD8+ T cells [#/µL] - N0 = Y(11) ;% resting non-SIV-specific CD8+ T cells [#/µL] - Na = Y(12)+Y(13) ;% active non-SIV-specific T cells [#/µL] + V = Y(01) ;% SIV virions [#/uL] + S0 = Y(02) ;% resting specific CD8+ T cells [#/uL] + Sa = sum(Y(03:10)) ;% active specific CD8+ T cells [#/uL] + N0 = Y(11) ;% resting non-SIV-specific CD8+ T cells [#/uL] + Na = Y(12)+Y(13) ;% active non-SIV-specific T cells [#/uL] C = Y(15) ;% N803 plasma concentration [pM] F = terms(Y') ;% modifiers for 'gS0','gN0','p0','aS0','aN0' p = F(13) ;% S0/N0 proliferation rate [/d] aS = F(14) ;% S0 activation rate [/d] aN = F(15) ;% N0 activation rate [/d] - gS = F(16) ;% Sa killing rate [#/µL-d] - gN = F(17) ;% Na killing rate [#/µL-d] + gS = F(16) ;% Sa killing rate [#/uL-d] + gN = F(17) ;% Na killing rate [#/uL-d] J = zeros(17) ; @@ -400,8 +399,8 @@ function F = terms(Y) - V = Y(:,01) ;% SIV virions [#/µL] - SN = sum(Y(:,02:13),2) ;% CD8+ T cells [#/µL] + V = Y(:,01) ;% SIV virions [#/uL] + SN = sum(Y(:,02:13),2) ;% CD8+ T cells [#/uL] C = Y(:,15) ;% N803 plasma concentration [pM] R = Y(:,17) ;% immune regulation [] @@ -427,8 +426,8 @@ F(:,13) = p0 * F(:,10).* F(:,02).* F(:,05) ;% S0/N0 proliferation [/d] F(:,14) = aS0* F(:,11).* F(:,03).* F(:,06) ;% S0 activation rate [/d] F(:,15) = aN0* (F(:,12)+F(:,04)).* F(:,07) ;% N0 activation rate [/d] - F(:,16) = gS0* F(:,08) ;% Sa killing rate [#/µL-d] - F(:,17) = gN0* F(:,09) ;% Na killing rate [#/µL-d] + F(:,16) = gS0* F(:,08) ;% Sa killing rate [#/uL-d] + F(:,17) = gN0* F(:,09) ;% Na killing rate [#/uL-d] F(:,18) = sS* F(:,11).* F(:,03) ;% regulation generation by aS [/d] F(:,19) = sN* (F(:,12) + F(:,04)) ;% regulation generation by aN [/d] diff --git a/N803_plotter_2_parameters.m b/N803_plotter_2_parameters.m new file mode 100644 index 0000000..09b8e09 --- /dev/null +++ b/N803_plotter_2_parameters.m @@ -0,0 +1,114 @@ +% N803_parameters.m - makes violin plots of parameter distributions + +close all +clear all %#ok +addpath Analyzer +addpath Data + +load('N803_MCMC','C') + +%% Model Parameters ======================================================= + +% Retrieving parameter boundaries ----------------------------------------- +warning off ; load('N803_shared_Results','Results') ; warning on +% q vS vN h pS pN d dA C50 aS1 aN1 dR p2 gN2 +id1 = [ 7 9 10 13 15 16 17 18 23 25 26 27 29 30 ] ; +id2 = [ 1 4 5 10 12 13 14 15 16 18 19 22 25 24 ] ; +Bounds = NaN*ones(2,27) ; +Bounds(:,id2) = Results(1).setup.GuessBOUNDS(:,id1) ; +Bounds(:,4:5) = Bounds(:,4:5)/1000 ;% converting vS,vN to [#/uL] + +% Plotting shared parameters ---------------------------------------------- +PARS = C(1).ParSample(:,[1:15,20:31]) ;% parameter sample +Size = [6.5 3] ;% figure size +figName = 'Fig S1 Parameters' ; + +f = figure('Name', figName, ... + 'Color', [1 1 1], ... + 'PaperUnits', 'Inches', ... + 'PaperSize', Size, ... + 'Units', 'Inches', ... + 'Position', [0 0 Size]) ; + +Names = {'q','g_S','g_N','V_{50,S}','V_{50,N}','a_S','a_N','m_S','m_N',... + 'H','p','p_S','p_N','d','d_A','C_{50}','\rho','\alpha_S','\alpha_N',... + 's_S','s_N','d_R','\lambda_S','\lambda_N','\phi','\zeta_S','\zeta_N'} ; + +parplotter(PARS, ... + 'AxesLims', [0 28 -3 4],... + 'AxesPos', [.075 .2 .9 .75],... + 'Bounds', Bounds,... + 'Factors', [1 1e3 1e3 ones(1,24)],...converting gS0,gN0 to [nL/#-d] + 'figHand', f,... + 'logScale', 1,... + 'Names', Names,... + 'ViolColor', [.5 .5 .5]) ; + +set(gca,'FontSize',10,'FontName','Arial','LineWidth',1,'XTickLabelRotation',315) +set(gcf,'InvertHardcopy','off') +print(figName,'-dtiff','-r300') + +%% Sampled Parameter Convolutions ========================================= +PARS = C(1).ParSample(:,1:31) ;% parameter sample + +% Calculating sampled parameter convolutions ------------------------------ +YIC1 = C(1).ParSample(:,32:end) ;% initial conditions +YIC2 = C(2).ParSample(:,32:end) ;% initial conditions +SN1 = sum( YIC1(:,2:13),2 ) ;% initial: S+N (co 1) +SN2 = sum( YIC2(:,2:13),2 ) ;% initial: S+N (co 2) +R1 = YIC1(:,end) ;% initial R (co 1) +R2 = YIC2(:,end) ;% initial R (co 2) +mS = PARS(:,08) ;% Sa reversion rate constant [/d] +mN = PARS(:,09) ;% Na reversion rate constant [/d] +h = PARS(:,10) ;% 50% S+N proliferation saturation [#/uL] +p0 = PARS(:,11) ;% S0/N0 proliferation rate constant [/d] +pS = PARS(:,12) ;% Sa proliferation rate constant [/d] +pN = PARS(:,13) ;% Na proliferation rate constant [/d] +d = PARS(:,14) ;% S0/N0 death rate constant [/d] +dA = PARS(:,15) ;% Sa/Na death rate constant [/d] +p1 = PARS(:,21) ;% S0/N0 proliferation stimulation factor [] +p2 = PARS(:,29) ;% S0/N0 proliferation regulation factor [] + +psi = PARS(:,03)./PARS(:,02) ;% gN0/gS0 +sig = PARS(:,25)./PARS(:,24) ;% sN/sS +fS1 = sum( YIC1(:,2:10),2 ) ./ SN1 ;% initial: S/(S+N) (co 1) +fS2 = sum( YIC2(:,2:10),2 ) ./ SN2 ;% initial: S/(S+N) (co 2) +fSn = (fS2-fS1)./ (.3-fS1) ;% normalized initial S/(S+N) (co 2) +US = 2*(2*pS./(pS+dA)).^7 ; +UN = 2* 2*pN./(pN+dA) ; +mSn = mS.*(US-1)./dA ;% normalized Sa reversion rate constant [] +mNn = mN.*(UN-1)./dA ;% normalized Na reversion rate constant [] +p0_max1 = d.*(1+p2.*R1).*(h+SN1)./h ; +p0_max2 = d.*(1+p2.*R2).*(h+SN2)./h ; +p0_max = min( [p0_max1 p0_max2],[],2 ) ; +p0n = p0./ p0_max ;% nomalized S0/N0 proliferation rate constant [/d] +pm = p1.*p0 ;% S0/N0 maximum proliferation rate [] + +% Plotting sampled parameter convolutions --------------------------------- +PARCONS = [ fS1 fSn psi mSn mNn p0n pm sig ] ; +id1 = [ 05 06 08 11 12 14 24 28 ] ; +Bounds = Results(1).setup.GuessBOUNDS(:,id1) ; +Size = [3 3*5/7] ;% figure size +figName = 'Fig S2 Convolutions' ;% figure name + +Names = {'f_S','\xi_S','g_N/g_S','\mu_S','\mu_N','\psi','p\cdot\rho','s_N/s_S'} ; + +f = figure('Name', figName, ... + 'Color', [1 1 1], ... + 'PaperUnits', 'Inches', ... + 'PaperSize', Size, ... + 'Units', 'Inches', ... + 'Position', [0 0 Size]) ; + +parplotter(PARCONS,... + 'AxesLims', [0 9 -2 2],... + 'AxesPos', [.2 .25 .75 .7],... + 'Bounds', Bounds,... + 'figHand', f,... + 'logScale', 1,... + 'Names', Names,... + 'ViolColor', [.5 .5 .5]) ; + +set(gca,'FontSize',10,'FontName','Arial','LineWidth',1,'XTickLabelRotation',315) +set(gcf,'InvertHardcopy','off') +print(figName,'-dtiff','-r300') diff --git a/N803_shared.m b/N803_shared.m index de3be4f..378bffe 100644 --- a/N803_shared.m +++ b/N803_shared.m @@ -1,19 +1,19 @@ %% N803_shared.m - solves model for 2 cohorts with shared parameters % % /--------------------------------------------------------------\ -% | Date: 12/25/2022 | +% | Date: 01/08/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 [#/ç„¡] -% S0 = resting SIV-specific CD8+ T cells [#/ç„¡] -% Sa = active SIV-specific CD8+ T cells [#/ç„¡] -% N0 = resting non-SIV-specific CD8+ T cells [#/ç„¡] -% Na = active non-SIV-specific CD8+ T cells [#/ç„¡] +% Nomenclature: V = SIV virions [#/uL] +% T8 = total CD8+ T cells [#/uL] +% S0 = resting SIV-specific CD8+ T cells [#/uL] +% Sa = active SIV-specific CD8+ T cells [#/uL] +% N0 = resting non-SIV-specific CD8+ T cells [#/uL] +% Na = active non-SIV-specific CD8+ T cells [#/uL] % X = N803 at absorption site [pmol/kg] % C = N803 plasma concentration [pM] % R = regulation [] (dimensionless quantity) @@ -75,8 +75,8 @@ % 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) -SNi(1) = AllPars(03) ;% S+N initial value [#/ç„¡] (cohort 1) -SNi(2) = AllPars(04) ;% S+N initial value [#/ç„¡] (cohort 2) +SNi(1) = AllPars(03) ;% S+N initial value [#/uL] (cohort 1) +SNi(2) = AllPars(04) ;% S+N initial value [#/uL] (cohort 2) fS(1) = AllPars(05) ;% initial frequency: S/(S+N) (cohort 1) fSn = AllPars(06) ;% normalized S/(S+N) (co 2) fS(2) = fS(1)+(0.3-fS(1))*fSn ;% initial frequency: S/(S+N) (cohort 2) @@ -88,7 +88,7 @@ mSn = AllPars(11) ;% normalized Sa reversion rate constant [] mNn = AllPars(12) ;% normalized Na reversion rate constant [] -SN50 = AllPars(13) ;% 50% S+N proliferation saturation [#/ç„¡] +SN50 = AllPars(13) ;% 50% S+N proliferation saturation [#/uL] p0n = AllPars(14) ;% nomalized S0/N0 proliferation rate constant [/d] pS = AllPars(15) ;% Sa proliferation rate constant [/d] pN = AllPars(16) ;% Na proliferation rate constant [/d] @@ -112,9 +112,9 @@ %% ------------------------------------------------------------------------ % Calculate some initial conditions & parameters -------------------------- -Vi = 10.^(Vi - 3) ;% converting V initial value [#/ç„¡] -V50S = V50S/1000 ;% 50% viral stimulation saturation for S [#/ç„¡] -V50N = V50N/1000 ;% 50% viral stimulation saturation for N [#/ç„¡] +Vi = 10.^(Vi - 3) ;% converting V initial value [#/uL] +V50S = V50S/1000 ;% 50% viral stimulation saturation for S [#/uL] +V50N = V50N/1000 ;% 50% viral stimulation saturation for N [#/uL] YS = Vi ./ (Vi+ V50S) ;% initial V/(V50S+V) [cohort 1,2] YN = Vi ./ (Vi+ V50N) ;% initial V/(V50N+V) [cohort 1,2] @@ -166,8 +166,8 @@ beta = psi * NA ./ (1+gN2*Ri) ; z = beta(1) - beta(2) ; a = z*R ; b = SA(1)*R - SA(2) + z*(1+R) ; c = SA(1) - SA(2) + z ; gS2 = ( -b + sqrt(b^2 - 4*a*c) ) / (2*a) ;% S killing regulation factor [] -gS0 = q / ( SA(1)/(1+gS2) + beta(1) ) ;% Sa killing rate constant [ç„¡/#-d] -gN0 = psi * gS0 ;% Na killing rate constant [ç„¡/#-d] +gS0 = q / ( SA(1)/(1+gS2) + beta(1) ) ;% Sa killing rate constant [uL/#-d] +gN0 = psi * gS0 ;% Na killing rate constant [uL/#-d] %% Do for each cohort (NOT indenting loop) ================================ for c = 1:2 @@ -187,17 +187,17 @@ % Prepare parameter and initial value vectors and call 'N803_model_2' ----- Pars(01) = q ;% V growth rate (if S+N were absent) [/d] -Pars(02) = gS0 ;% Sa killing rate constant [ç„¡/#-d] -Pars(03) = gN0 ;% Na killing rate constant [ç„¡/#-d] +Pars(02) = gS0 ;% Sa killing rate constant [uL/#-d] +Pars(03) = gN0 ;% Na killing rate constant [uL/#-d] -Pars(04) = V50S ;% 50% viral stimulation saturation for S [#/ç„¡] -Pars(05) = V50N ;% 50% viral stimulation saturation for N [#/ç„¡] +Pars(04) = V50S ;% 50% viral stimulation saturation for S [#/uL] +Pars(05) = V50N ;% 50% viral stimulation saturation for N [#/uL] Pars(06) = aS0 ;% S0 activation rate constant [/d] Pars(07) = aN0 ;% N0 activation rate constant [/d] Pars(08) = mS ;% Sa reversion rate constant [/d] Pars(09) = mN ;% Na reversion rate constant [/d] -Pars(10) = SN50 ;% 50% S+N proliferation saturation [#/ç„¡] +Pars(10) = SN50 ;% 50% S+N proliferation saturation [#/uL] Pars(11) = p0 ;% S0/N0 proliferation rate constant [/d] Pars(12) = pS ;% Sa proliferation rate constant [/d] Pars(13) = pN ;% Na proliferation rate constant [/d] @@ -252,4 +252,4 @@ Y_OUT = [ Y_Cohort{1} , Y_Cohort{2} ] ; PARS = [ P_Cohort{1} ; P_Cohort{2} ] ; -end +end \ No newline at end of file diff --git a/N803_single.m b/N803_single.m index b42ec48..ce742b0 100644 --- a/N803_single.m +++ b/N803_single.m @@ -1,19 +1,19 @@ %% N803_single.m - solves model for 1 cohort % % /--------------------------------------------------------------\ -% | Date: 12/25/2022 | +% | Date: 01/08/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 [#/ç„¡] -% S0 = resting SIV-specific CD8+ T cells [#/ç„¡] -% Sa = active SIV-specific CD8+ T cells [#/ç„¡] -% N0 = resting non-SIV-specific CD8+ T cells [#/ç„¡] -% Na = active non-SIV-specific CD8+ T cells [#/ç„¡] +% Nomenclature: V = SIV virions [#/uL] +% T8 = total CD8+ T cells [#/uL] +% S0 = resting SIV-specific CD8+ T cells [#/uL] +% Sa = active SIV-specific CD8+ T cells [#/uL] +% N0 = resting non-SIV-specific CD8+ T cells [#/uL] +% Na = active non-SIV-specific CD8+ T cells [#/uL] % X = N803 at absorption site [pmol/kg] % C = N803 plasma concentration [pM] % R = regulation [] (dimensionless quantity) @@ -59,7 +59,7 @@ % Rename inputed parameters ----------------------------------------------- Vi = AllPars(01) ;% V initial value [log(#/mL)] -SNi = AllPars(02) ;% S+N initial value [#/ç„¡] +SNi = AllPars(02) ;% S+N initial value [#/uL] fS = AllPars(03) ;% initial frequency: S/(S+N) fSA = AllPars(04) ;% initial frequency: Sa/S q = AllPars(05) ;% V growth rate (if S+N were absent) [/d] @@ -68,7 +68,7 @@ V50N = AllPars(08) ;% 50% viral stimulation saturation for N [#/mL] mSn = AllPars(09) ;% normalized Sa reversion rate constant [] mNn = AllPars(10) ;% normalized Na reversion rate constant [] -SN50 = AllPars(11) ;% 50% S+N proliferation saturation [#/ç„¡] +SN50 = AllPars(11) ;% 50% S+N proliferation saturation [#/uL] pS = AllPars(12) ;% Sa proliferation rate constant [/d] pN = AllPars(13) ;% Na proliferation rate constant [/d] d = AllPars(14) ;% S0/N0 death rate constant [/d] @@ -92,9 +92,9 @@ %% ------------------------------------------------------------------------ % Calculate some initial conditions & parameters -------------------------- -Vi = 10^(Vi-3) ;% V initial value [#/ç„¡] -V50S = V50S/1000 ;% 50% viral stimulation saturation for S [#/ç„¡] -V50N = V50N/1000 ;% 50% viral stimulation saturation for N [#/ç„¡] +Vi = 10^(Vi-3) ;% V initial value [#/uL] +V50S = V50S/1000 ;% 50% viral stimulation saturation for S [#/uL] +V50N = V50N/1000 ;% 50% viral stimulation saturation for N [#/uL] % restrict mS and mN such that initial activation aS and aN are positive US = 2*(2*pS/(pS+dA))^7 ; @@ -146,8 +146,8 @@ p0 = pi * (SN50+SNi)/SN50 * (1+p2) ;% S0/N0 proliferation rate con [/d] p1 = pm/p0 ;% S0/N0 proliferation stimulation factor -gS0 = q / ( SA/(1+gS2) + psi*NA/(1+gN2) ) ;% Sa killing rate con [ç„¡/#-d] -gN0 = psi*gS0 ;% Na killing rate constant [ç„¡/#-d] +gS0 = q / ( SA/(1+gS2) + psi*NA/(1+gN2) ) ;% Sa killing rate con [uL/#-d] +gN0 = psi*gS0 ;% Na killing rate constant [uL/#-d] sS = dR / ( Vi/(V50S+Vi) + sig*Vi/(V50N+Vi) ) ;% regulation gen by S act sN = sig*sS ;% regulation generation by N activation @@ -156,17 +156,17 @@ % Prepare parameter and initial value vectors and call 'N803_model_2' ----- Pars(01) = q ;% V growth rate (if S+N were absent) [/d] -Pars(02) = gS0 ;% Sa killing rate constant [ç„¡/#-d] -Pars(03) = gN0 ;% Na killing rate constant [ç„¡/#-d] +Pars(02) = gS0 ;% Sa killing rate constant [uL/#-d] +Pars(03) = gN0 ;% Na killing rate constant [uL/#-d] -Pars(04) = V50S ;% 50% viral stimulation saturation for S [#/ç„¡] -Pars(05) = V50N ;% 50% viral stimulation saturation for N [#/ç„¡] +Pars(04) = V50S ;% 50% viral stimulation saturation for S [#/uL] +Pars(05) = V50N ;% 50% viral stimulation saturation for N [#/uL] Pars(06) = aS0 ;% S0 activation rate constant [/d] Pars(07) = aN0 ;% N0 activation rate constant [/d] Pars(08) = mS ;% Sa reversion rate constant [/d] Pars(09) = mN ;% Na reversion rate constant [/d] -Pars(10) = SN50 ;% 50% S+N proliferation saturation [#/ç„¡] +Pars(10) = SN50 ;% 50% S+N proliferation saturation [#/uL] Pars(11) = p0 ;% S0/N0 proliferation rate constant [/d] Pars(12) = pS ;% Sa proliferation rate constant [/d] Pars(13) = pN ;% Na proliferation rate constant [/d] @@ -214,4 +214,4 @@ Pars = [ Pars Yic ] ;% adding initial conditions to parameter output -end +end \ No newline at end of file diff --git a/main_shared_calibration_2022_08_18.m b/main_shared_calibration_2022_08_18.m index ca3409e..526c7ee 100644 --- a/main_shared_calibration_2022_08_18.m +++ b/main_shared_calibration_2022_08_18.m @@ -37,8 +37,8 @@ % GuessBOUNDS = array for boundaries of initial parameter guesses ----- GuessBOUNDS(:,01) = [1 1]*mean(cat(1,IVs{1}{:})) ;% V iv [log(#/mL)] (co 1) GuessBOUNDS(:,02) = [1 1]*mean(cat(1,IVs{3}{:})) ;% V iv [log(#/mL)] (co 2) -GuessBOUNDS(:,03) = [1 1]*mean(cat(1,IVs{2}{:})) ;% S+N iv [#/µL] (co 1) -GuessBOUNDS(:,04) = [1 1]*mean(cat(1,IVs{4}{:})) ;% S+N iv [#/µL] (co 2) +GuessBOUNDS(:,03) = [1 1]*mean(cat(1,IVs{2}{:})) ;% S+N iv [#/uL] (co 1) +GuessBOUNDS(:,04) = [1 1]*mean(cat(1,IVs{4}{:})) ;% S+N iv [#/uL] (co 2) GuessBOUNDS(:,05) = [.03 .3] ;% initial frequency: S/(S+N) (co 1) GuessBOUNDS(:,06) = [.01 1] ;% normalized S/(S+N) (co 2) @@ -49,7 +49,7 @@ GuessBOUNDS(:,11) = [1 10] ;% normalized Sa reversion rate constant [] GuessBOUNDS(:,12) = [1 10] ;% normalized Na reversion rate constant [] -GuessBOUNDS(:,13) = [300 3e3] ;% 50% S+N proliferation saturation [#/µL] +GuessBOUNDS(:,13) = [300 3e3] ;% 50% S+N proliferation saturation [#/uL] GuessBOUNDS(:,14) = [0.1 1] ;% S0/N0 normalized prolif rate constant [/d] GuessBOUNDS(:,15) = [1 4] ;% Sa proliferation rate constant [/d] GuessBOUNDS(:,16) = [0.5 2] ;% Na proliferation rate constant [/d] diff --git a/main_singles_calibration_2022_09_25.m b/main_singles_calibration_2022_09_25.m index 2fdc22a..e74bb6e 100644 --- a/main_singles_calibration_2022_09_25.m +++ b/main_singles_calibration_2022_09_25.m @@ -33,7 +33,7 @@ GuessBOUNDS(:,08) = [1e4 1e6] ;% 50% viral stim saturation for N [#/mL] GuessBOUNDS(:,09) = [1 10] ;% normalized Sa reversion rate constant [] GuessBOUNDS(:,10) = [1 10] ;% normalized Na reversion rate constant [] -GuessBOUNDS(:,11) = [300 3e3] ;% 50% S+N proliferation saturation [#/µL] +GuessBOUNDS(:,11) = [300 3e3] ;% 50% S+N proliferation saturation [#/uL] GuessBOUNDS(:,12) = [1 4] ;% Sa proliferation rate constant [/d] GuessBOUNDS(:,13) = [0.5 2] ;% Na proliferation rate constant [/d] GuessBOUNDS(:,14) = [.01 .05] ;% S0/N0 death rate constant [/d] @@ -77,7 +77,7 @@ % GuessBOUNDS = array for boundaries of initial parameter guesses --------- GuessBOUNDS(:,01) = [1 1]*mean(cat(1,IVs{1}{:})) ;% V iv [log(#/mL)] -GuessBOUNDS(:,02) = [1 1]*mean(cat(1,IVs{2}{:})) ;% S0+Sa iv [#/µL] +GuessBOUNDS(:,02) = [1 1]*mean(cat(1,IVs{2}{:})) ;% S0+Sa iv [#/uL] % FitBOUNDS = array for boundaries of parameter fitting ------------------- FitBOUNDS = GuessBOUNDS ;