I have some errors in executing the following code, can anyone help me please?

% Define PK model parameters

infusion_rate = 200; % ug/kg/min

CL1 = 0.066; % L/min

CL2 = 0.6; % L/min

CL3 = 0.35; % L/min

V1 = 3; % L

V2 = 20; % L

V3 = 60; % L

% Define PD model parameters

ke0 = 0.19; % 1/min

ke1 = 0.07; % 1/min

ke2 = 0.04; % 1/min

EC50 = 2.2; % ug/mL

% Define time span and initial drug concentrations

tspan = [0 30]; % minutes

Cp0 = [0 0 0]; % ug/mL

Ce0 = 0; % ug/mL

% Define PK-PD model as an anonymous function

propofol_model = @(t, y) [ (infusion_rate - CL1 * y(1) - CL2 *(y(1) - y(2))/V1 - CL3*(y(1) - y(3))/V1)/V1; ...

(CL2*(y(1) - y(2))/V2); ...

(CL3*(y(1) - y(3))/V3)]; % PK model

pd_model = @(t, y) [ -ke0 * y(1) + ke1 * y(4); ... % Effect site 1

-ke1 * y(2) + ke2 * y(5); ... % Effect site 2

-ke2*y(3); ... % Effect site 3

(y(1) + y(2) + y(3))/EC50 - y(7)]; % PD model

% Solve the ODEs using ode45 solver

options = odeset('RelTol',1e-6,'AbsTol',1e-6); % Set solver options

[t, y] = ode45(@(t, y) [propofol_model(t, y(1:3)); pd_model(t, y)], tspan, [Cp0 Ce0], options);

% Plot PK and PD results

figure;

subplot(2,1,1);

plot(t, y(:,1), 'r', t, y(:,2), 'g', t, y(:,3), 'b');

xlabel('Time (min)');

ylabel('Drug concentration (ug/mL)');

legend('Central compartment', 'Peripheral compartment 1', 'Peripheral compartment 2');

title('Propofol PK Model');

subplot(2,1,2);

plot(t, y(:,7));

xlabel('Time (min)');

ylabel('Effect site concentration (ug/mL)');

title('Propofol PD Model');

Answer 1

[Alan Stevens]

You have no initial guesses for y(5) and y(7). And what about y(6)?

Answer 2

[Torsten]

The vector of initial values

[Cp0 Ce0]

must have 7 elements, not 4.