Skip to content

Commit

Permalink
Add files via upload
Browse files Browse the repository at this point in the history
  • Loading branch information
@cody2
cody2 authored Jul 8, 2022
1 parent 7461fcb commit 5791e28
Showing 1 changed file with 93 additions and 71 deletions.
164 changes: 93 additions & 71 deletions N803_model_2.m
View file
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
%% N803_model_2.m - solves ODE model for N-803 treatment of SIV
%
% /--------------------------------------------------------------\
% | Date: 06/29/2022 |
% | Date: 07/08/2022 |
% | Author: Jonathan Cody |
% | Affiliation: Purdue University |
% | Weldon School of Biomedical Engineering |
Expand Down Expand Up @@ -165,67 +165,86 @@

%% ------------------------------------------------------------------------
% Add columns to 'Y_OUT' (if requested) -----------------------------------
% NOTE: ERC1 = effective first-order growth rate constant

if isempty(fullOut) ; return ; end

Y_OUT = [ Y_OUT, Y(:,2:end), zeros(Pieces(end),33) ] ;

F = terms(Y) ;% see 'terms' function
gE = gE0*F(:,08) ;% specific killing rate [#/無-d]
gB = gB0*F(:,09) ;% bystander killing rate [#/無-d]
p = p0*F(:,13) ;% E0/B0 proliferation rate [/d]
aE = aE0*F(:,14) ;% specific activation rate [/d]
aB = aB0*F(:,15) ;% bystander activation rate [/d]

Y_OUT(:,19) = F(:,01) ;% N803 effect []
Y_OUT = [ Y_OUT, Y, zeros(Pieces(end),46) ] ;

V = Y(:,01) ;% SIV virions [#/無]
E0 = Y(:,02) ;% resting specific CD8+ T cells [#/無]
B0 = Y(:,11) ;% resting bystander CD8+ T cells [#/無]
Ea = sum( Y(:,03:10) ,2) ;% active specific CD8+ T cells [#/無]
Ba = sum( Y(:,12:13) ,2) ;% active bystander T cells [#/無]
Y_OUT(:,20) = Ea ;% Ea
Y_OUT(:,21) = Ba ;% Ba
Y_OUT(:,22) = (E0+Ea) ;% E
Y_OUT(:,23) = (B0+Ba) ;% B
Y_OUT(:,24) = (E0+Ea)./(E0+B0+Ea+Ba) ;% E frequency in (E+B)
Y_OUT(:,25) = (Ea)./(E0+Ea) ;% activated frequency in E
Y_OUT(:,26) = (Ba)./(B0+Ba) ;% activated frequency in B

Eaf = Y(:,10) ;% last generation of Ea
Baf = Y(:,13) ;% last generation of Ba
Y_OUT(:,27) = p - d ;% E0/B0 expansion ERC1
Y_OUT(:,28) = mE*Eaf./E0 - aE ;% E0 activation ERC1
Y_OUT(:,29) = mB*Baf./B0 - aB ;% B0 activation ERC1
Y_OUT(:,30) = log10( F(:,02) ) ;% C multiplier for p* (log)
Y_OUT(:,31) = log10( F(:,05) ) ;% R multiplier for p* (log)
Y_OUT(:,32) = log10( F(:,10) ) ;% density-dep multiplier for p* (log)

Eap = sum( Y(:,03:09) ,2) ;% Ea that is proliferating
Bap = Y(:,12) ;% Ba that is proliferating
Y_OUT(:,33) = (pE*Eap - dA*Ea) ./ Ea ;% EA expansion ERC1
Y_OUT(:,34) = (pB*Bap - dA*Ba) ./ Ba ;% BA expansion ERC1
Y_OUT(:,35) = (aE.*E0 - mE*Eaf) ./ Ea ;% EA activation ERC1
Y_OUT(:,36) = (aB.*B0 - mB*Baf) ./ Ba ;% BA activation ERC1
Y_OUT(:,37) = log10( F(:,03) ) ;% C multiplier for aE* (log)
Y_OUT(:,38) = log10( ( F(:,04) + F(:,12) )./F(:,12) ) ;% same for aB* (log)
Y_OUT(:,39) = log10( F(:,06) ) ;% R multiplier for aE* (log)
Y_OUT(:,40) = log10( F(:,07) ) ;% R multiplier for aB* (log)
Y_OUT(:,41) = log10( F(:,11) ) ;% V multiplier for aE* (log)
Y_OUT(:,42) = log10( F(:,12) ) ;% V multiplier for aB* (log)

Y_OUT(:,43) = q - gE.*Ea - gB.*Ba ;% V ERC1
Y_OUT(:,44) = log10( gE.*Ea ) ;% E killing ERC1 (log)
Y_OUT(:,45) = log10( gB.*Ba ) ;% B killing ERC1 (log)
Y_OUT(:,46) = log10( F(:,08) ) ;% R multiplier for gE* (log)
Y_OUT(:,47) = log10( F(:,09) ) ;% R multiplier for gB* (log)
Y_OUT(:,48) = gB.*Ba ./ (gE.*Ea + gB.*Ba) ;% B killing fraction

RgenE = sE* F(:,11).* F(:,03) ;% R gen by aE
RgenB = sB*( F(:,12) + F(:,04) ) ;% R gen by aB
Y_OUT(:,49) = log10( RgenE ) ;% R gen by aE (log)
Y_OUT(:,50) = log10( RgenB ) ;% R gen by aB (log)
Y_OUT(:,51) = RgenB ./ (RgenE + RgenB) ;% R gen by aB fraction
Eaf = Y(:,10) ;% last generation of Ea
Baf = Y(:,13) ;% last generation of Ba
Eap = Ea - Eaf ;% Ea that is proliferating
Bap = Y(:,12) ;% Ba that is proliferating

F = terms(Y) ;% see 'terms' function
p = F(:,13) ;% E0/B0 proliferation rate [/d]
aE = F(:,14) ;% E0 activation rate [/d]
aB = F(:,15) ;% B0 activation rate [/d]
gE = F(:,16) ;% Ea killing rate [#/無-d]
gB = F(:,17) ;% Ba killing rate [#/無-d]

Y_OUT(:,20) = F(:,01) ;% N803 effect []

Y_OUT(:,21) = (E0+B0+Ea+Ba) ;% CD8+ T cells [#/無]
Y_OUT(:,22) = Ea ;% Ea
Y_OUT(:,23) = Ba ;% Ba
Y_OUT(:,24) = (E0+Ea) ;% E
Y_OUT(:,25) = (B0+Ba) ;% B
Y_OUT(:,26) = (E0+Ea)./(E0+B0+Ea+Ba) ;% E frequency in (E+B)
Y_OUT(:,27) = (Ea)./(E0+Ea) ;% activated frequency in E
Y_OUT(:,28) = (Ba)./(B0+Ba) ;% activated frequency in B
Y_OUT(:,29) = log10( Ea ) ;% Ea (log)
Y_OUT(:,30) = log10( Ba ) ;% Ba (log)

Y_OUT(:,31) = p.*E0 ;% E0 proliferation term [#/無/d]
Y_OUT(:,32) = p.*B0 ;% B0 proliferation term [#/無/d]
Y_OUT(:,33) = d.*E0 ;% E0 death term [#/無/d]
Y_OUT(:,34) = d.*B0 ;% B0 death term [#/無/d]
Y_OUT(:,35) = p ;% E0/B0 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) = aE.*E0 ;% E0 activation term [#/無/d]
Y_OUT(:,40) = aB.*B0 ;% B0 activation term [#/無/d]
Y_OUT(:,41) = pE*Eap ;% Ea proliferation term [#/無/d]
Y_OUT(:,42) = pB*Bap ;% Ba proliferation term [#/無/d]
Y_OUT(:,43) = dA*Ea ;% Ea death term [#/無/d]
Y_OUT(:,44) = dA*Ba ;% Ba death term [#/無/d]
Y_OUT(:,45) = mE*Eaf ;% Ea reversion term [#/無/d]
Y_OUT(:,46) = mB*Baf ;% Ba reversion term [#/無/d]
Y_OUT(:,47) = aE ;% E0 activation rate [/d]
Y_OUT(:,48) = aB ;% B0 activation rate [/d]
Y_OUT(:,49) = F(:,03) ;% C multiplier for aE*
Y_OUT(:,50) = ( F(:,04) + F(:,12) )./F(:,12) ;% C multiplier for aE*
Y_OUT(:,51) = F(:,06) ;% R multiplier for aE*
Y_OUT(:,52) = F(:,07) ;% R multiplier for aB*
Y_OUT(:,53) = F(:,11) ;% V multiplier for aE*
Y_OUT(:,54) = F(:,12) ;% V multiplier for aB*

Y_OUT(:,55) = q*V ;% V growth term [#/無/d]
Y_OUT(:,56) = gE.*Ea.*V ;% E killing term [#/無/d]
Y_OUT(:,57) = gB.*Ba.*V ;% B killing term [#/無/d]
Y_OUT(:,58) = gE.*Ea ;% E killing rate [/d]
Y_OUT(:,59) = gB.*Ba ;% B killing rate [/d]
Y_OUT(:,60) = F(:,08) ;% R multiplier for gE*
Y_OUT(:,61) = F(:,09) ;% R multiplier for gB*
Y_OUT(:,62) = gB.*Ba ./ (gE.*Ea + gB.*Ba) ;% B killing fraction

Y_OUT(:,63) = F(:,18) ;% R gen by aE
Y_OUT(:,64) = F(:,19) ;% R gen by aB
Y_OUT(:,65) = F(:,19) ./ (F(:,18) + F(:,19)) ;% R gen by aB fraction

% Converting most outputs to log-value
logID = [ 3:15 , 18:19 , 31:61 , 63:64 ] ;
Y_OUT(:,logID) = log10( Y_OUT(:,logID) ) ;

% Converting virus to #/mL
Y_OUT(:,3) = Y_OUT(:,3) + 3 ;

%% ========================================================================
% NESTED FUNCTIONS
Expand All @@ -247,11 +266,12 @@
R = Y(16:17) ;% immune regulation []
F = terms(Y') ;% modifiers for 'gE0','gB0','p0','aE0','aB0'

gE = gE0*F(08) ;% specific killing rate [#/無-d]
gB = gB0*F(09) ;% bystander killing rate [#/無-d]
p = p0*F(13) ;% E0/B0 proliferation rate [/d]
aE = aE0*F(14) ;% specific activation rate [/d]
aB = aB0*F(15) ;% bystander activation rate [/d]
p = F(13) ;% E0/B0 proliferation rate [/d]
aE = F(14) ;% E0 activation rate [/d]
aB = F(15) ;% B0 activation rate [/d]
gE = F(16) ;% Ea killing rate [#/無-d]
gB = F(17) ;% Ba killing rate [#/無-d]
s = F(18)+F(19) ;% regulation generation rate [/d]

dY = Y ;% dY/dt

Expand All @@ -266,7 +286,7 @@
dY(13) = 2*pB*Ba(1) - dA*Ba(2) - mB*Ba(2) ;% dB2/dt
dY(14) = - ka*X ;% dX/dt
dY(15) = ka*X/vd - ke*C ;% dC/dt
dY(16) = F(16) - dR*R(1) ;% dR1/dt
dY(16) = s - dR*R(1) ;% dR1/dt
dY(17) = dR*R(1) - dR*R(2) ;% dR2/dt

end
Expand All @@ -281,11 +301,11 @@
C = Y(15) ;% N803 plasma concentration [pM]
F = terms(Y') ;% modifiers for 'gE0','gB0','p0','aE0','aB0'

gE = gE0*F(08) ;% specific killing rate [#/無-d]
gB = gB0*F(09) ;% bystander killing rate [#/無-d]
p = p0*F(13) ;% E0/B0 proliferation rate [/d]
aE = aE0*F(14) ;% specific activation rate [/d]
aB = aB0*F(15) ;% bystander activation rate [/d]
p = F(13) ;% E0/B0 proliferation rate [/d]
aE = F(14) ;% E0 activation rate [/d]
aB = F(15) ;% B0 activation rate [/d]
gE = F(16) ;% Ea killing rate [#/無-d]
gB = F(17) ;% Ba killing rate [#/無-d]

J = zeros(17) ;

Expand Down Expand Up @@ -372,7 +392,7 @@
R = Y(:,17) ;% immune regulation []

% percolate vector for storing terms
F = zeros(size(Y,1),16) ;
F = zeros(size(Y,1),19) ;

F(:,01) = C./ (C50+C) ;% N803 effect []

Expand All @@ -390,12 +410,14 @@
F(:,11) = V./ (V50E+V) ;% Ea activation antigen dependence []
F(:,12) = V./ (V50B+V) ;% Ba activation antigen dependence []

F(:,13) = F(:,10).* F(:,02).* F(:,05) ;% E0/B0 prolif modifier []
F(:,14) = F(:,11).* F(:,03).* F(:,06) ;% Ea activation modifier []
F(:,15) = (F(:,12)+F(:,04)).* F(:,07) ;% Ba activation modifier []
F(:,13) = p0 * F(:,10).* F(:,02).* F(:,05) ;% E0/B0 proliferation [/d]
F(:,14) = aE0* F(:,11).* F(:,03).* F(:,06) ;% E0 activation rate [/d]
F(:,15) = aB0* (F(:,12)+F(:,04)).* F(:,07) ;% B0 activation rate [/d]
F(:,16) = gE0* F(:,08) ;% Ea killing rate [#/無-d]
F(:,17) = gB0* F(:,09) ;% Ba killing rate [#/無-d]

F(:,16) = sE* F(:,11).* F(:,03) ...
+ sB*( F(:,12) + F(:,04) ) ;% regulation gen [/d]
F(:,18) = sE* F(:,11).* F(:,03) ;% regulation generation by aE [/d]
F(:,19) = sB* (F(:,12) + F(:,04)) ;% regulation generation by aB [/d]
end

end

0 comments on commit 5791e28

Please sign in to comment.