Trying to write a 3 compartment model
5 次查看(过去 30 天)
显示 更早的评论
% set parameter values
Vd = 210; % volume of distribution litres = vol compartment 2
V1 = 10; % volume of gut compartment
V3 = 30;
Dose = 25.*1000; % milligrams* 1000 gives micrograms
C0 = Dose/V1;
k1 = 0.24; % per hour
k2 = 0.28; % per hour
vmax = .47;
km = .38;
a = .42;
tmax = 12; % time range over which to solve model (hours)
Cinit = [C0,0,0]; % initial value for compartment 1 and compartment 2
% write down the expressions for the fluxes
J1i = @(t,C) 0; % no continuous influx to compartment 1
J1e = @(t,C) k1.*V1.*C(1); % efflux from compartment 1
J2i = J1e;
J2e = @(t,C) k2.*V2.*C(2); % influx to compartment 2 = efflux from compartment 1
J3i = J2e;
J3e = @(t,C) k2.*V3.*C(3); % efflux from compartment 3
% solve the differential equation
sol = ode45(@ThreeCompartmentModel,[0 tmax],Cinit,[],J1i,J1e,J2i,J2e,J3i,J3e,V1,Vd,V3);% plot results
% plot results
tstep = 0.01;
tspan = 0:tstep:tmax;
y = deval(sol,tspan);
plot(tspan,y)
legend('Compartment 1','Compartment 2','Concentration 3')
ylabel('drug concentration')
xlabel('time (hours)')
I wanted to do this in SimBiology or something. And be able to do multiple doses. But when I look at the SimBiology, it is modeling data I think. Can I simply represent a system like this? I wanted to model this. https://php.radford.edu/~ejmt/deliveryBoy.php?paper=eJMT_v1n3p3
2 个评论
Praveen Kumar M
2018-11-19
Hi Ryan,
I am clinical pharamcologist. Have you found the answer? If yes, please share the solution with us. And if you have not found the answer, then we can combinely work towards the solution. You can contact me either through MathWorks or by gmail (my mail: praveenkumarpgiindia@gmail.com)
Praveen Kumar M
2018-11-19
In the above code, which you have shared, I have a query. How have you defined the ThreeCompartmentModel function or variable?
回答(2 个)
Achal Mahajan
2022-6-1
Hi Ryan,
Here is the implementation of the problem in SimBiology. The following implementation is motivated from the paper that you have shared. To understand it better, as Arthur have mentioned above, please use the examples/videos/webinars.
clear all;clc;close all;
model = sbiomodel('Q4_MATLAB');
comp1Obj = addcompartment(model,'Stomach');
comp1Obj.CapacityUnits = 'liter';
comp1Obj.Value = 1;
comp2Obj = addcompartment(model,'Small_Int');
comp2Obj.CapacityUnits = 'liter';
comp2Obj.Value = 1;
comp3Obj = addcompartment(model,'Central');
comp3Obj.CapacityUnits = 'liter';
comp3Obj.Value = 1;
% Let's define all the reactions in their respective compartments
r1 = addreaction(model,'Stomach.Alcohol -> Small_Int.Alcohol');
r1.ReactionRate = 'k1*Stomach.Alcohol/(1+a*Stomach.Alcohol^2)';
r2 = addreaction(model,'Small_Int.Alcohol -> Central.Alcohol');
k = r2.addkineticlaw('MassAction');
k.ParameterVariableNames = 'k2';
r3 = addreaction(model,'Central.Alcohol -> null');
k = r3.addkineticlaw('Henri-Michaelis-Menten');
k.ParameterVariableNames = {'Vmax','km'};
k.SpeciesVariableNames = 'Central.Alcohol';
% Initial concentration and specify units
for i = 1:length(model.Species)
model.Species(i).Units = 'gram/liter';
end
model.Species(1).InitialAmount = 0.455;% This is what is mentioned in the paper, different from what you have implemented
model.Species(2).InitialAmount = 0;
model.Species(3).InitialAmount = 0;
% Now we can add the parameters of the model like rate constants
addparameter(model,'Vmax','Value',0.470,'ValueUnits','gram/liter/hour');
addparameter(model,'km','Value',0.380,'ValueUnits','gram/liter');
addparameter(model,'k1','Value',5.55,'ValueUnits','1/hour');
addparameter(model,'k2','Value',7.05,'ValueUnits','1/hour');
addparameter(model,'a','Value',0.42,'ValueUnits','liter^2/gram^2');
% Let's run the simulation now
cs = getconfigset(model,'active');
cs.StopTime = 4;
cs.TimeUnits = 'hour';
% Run the simulation and plot the data
cs.CompileOptions.UnitConversion = true;
[t, simdata, names] = sbiosimulate(model);
%[t, simdata, names] = sbiosimulate(model,d1); %Uncomment this if you want to add a dose
%[t, simdata, names] = sbiosimulate(model,v2); % or v3 Uncomment this if you want to add a variant
plot(t,simdata)
legend(names,'Location','NorthEastOutside');
-Achal
0 个评论
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Scan Parameter Ranges 的更多信息
产品
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!