How to select which parameters will be optimized in a model?

I am trying to do a simple optimization problem, but keep running into an error. The basic idea is shown in the code below.

% Extract parameters that need to be optimized
params_to_optimize = [true, true, true];

% Model parameters
k1_initial = 0.1;  % Initial guess for k1
k2_initial = 0.05; % Initial guess for k3
k3_initial = 0.02; % Initial guess for k2

% Concatenate initial guesses
initial_guess = [k1_initial, k2_initial, k3_initial];

% Bounds for optimization
lb = [0, 0, 0]; % Lower bounds for k1 and k3
ub = [1, 1, 1]; % Upper bounds for k1 and k3


% Provided experimental data
experimental_data_A = [1, 0.8, 0.7];
experimental_data_B = [60, 50, 48];
experimental_data_C = [100, 90, 80];

% Time points corresponding to the experimental data
experimental_time_points = [2, 4, 6];

% Initial conditions
A0 = 100;   % Initial concentration of A
B0 = 150;   % Initial concentration of B
C0 = 0;     % Initial concentration of C

% Time span
tspan = linspace(0, 6, 100);



% Objective function for patternsearch
objective_function = @(params) norm(simulate_and_compare(params, tspan, [A0, B0, C0], experimental_time_points, [experimental_data_A; experimental_data_B; experimental_data_C],params_to_optimize), 'fro');

% Options for the pattern search optimization
options_optimization = optimoptions('patternsearch', 'Display', 'iter', 'PlotFcn', @psplotbestf);

% Perform pattern search optimization
fitted_parameters = patternsearch(objective_function, initial_guess(params_to_optimize), [], [], [], [], lb(params_to_optimize), ub(params_to_optimize), [], options_optimization);

% Simulate the ODE system with the optimized parameters
[t_fit, y_fit] = ode15s(@(t, y) odeSystem(t, y, fitted_parameters), tspan, [A0, B0, C0]);

% % Plot the results including experimental data and fitted model
% figure;
% plot(t_fit, y_fit(:, 1), 'LineWidth', 2, 'DisplayName', 'A (Fitted Model)');
% hold on;
% plot(t_fit, y_fit(:, 2), 'LineWidth', 2, 'DisplayName', 'B (Fitted Model)');
% plot(t_fit, y_fit(:, 3), 'LineWidth', 2, 'DisplayName', 'C (Fitted Model)');
% scatter(experimental_time_points, experimental_data_A, 100, 'o', 'DisplayName', 'A (Experimental)');
% scatter(experimental_time_points, experimental_data_B, 100, 's', 'DisplayName', 'B (Experimental)');
% scatter(experimental_time_points, experimental_data_C, 100, '^', 'DisplayName', 'C (Experimental)');
% xlabel('Time');
% ylabel('Concentration');
% legend('Location', 'Best');
% title('Simulation and Experimental Data with Fitted Model');

% Display the optimized parameters
disp('Optimized Parameters:');
disp(fitted_parameters);

% ODE system
function dydt = odeSystem(t, y, params)
    k1 = params(1);
    k2 = params(2);
    k3 = params(3);
     

    % ODEs
    dydt = [
        -k1 * y(1) * y(2) + k3 * y(3);   % ODE for A
        -k1 * y(1) * y(2) + k3 * y(3);   % ODE for B
        k1 * y(1) * y(2) - k3 * y(3) - k2 * y(3)  % ODE for C
    ];
end

% Function to simulate the ODE system and compare with experimental data
function residuals = simulate_and_compare(params, tspan, initial_conditions, exp_time_points, exp_data, params_to_optimize)
    % Extract parameters to be optimized
    params_optimized = params(params_to_optimize);

    % Set the fixed parameters to their initial guesses
    initial_guess = initial_conditions;
    initial_guess(params_to_optimize) = params_optimized;

    % Simulate the ODE system
    [~, y_sim] = ode15s(@(t, y) odeSystem(t, y, initial_guess), tspan, initial_conditions);

    % Interpolate simulated data at experimental time points
    y_sim_interp = interp1(tspan, y_sim, exp_time_points);

    % Flatten the arrays for comparison
    y_sim_flat = y_sim_interp(:);
    exp_data_flat = exp_data(:);

    % Calculate residuals
    residuals = y_sim_flat - exp_data_flat;
end

In this case, I am asking it to optimize all three parameters in the model based on

params_to_optimize = [true, true, true];

However, if I change the params_to_optimize to other settings (eg: true, false, true), it returns an error.

The logical indices contain a true value outside of the array bounds.

Error in simulate_and_compare 
    params_optimized = params(params_to_optimize);

My goal is to be able to intially determine which of the parameters will be optimized (not necessarily all of them) and run the simulations accordingly. How can I do that?

Thanks!

Answer 1

[Sulaymon Eshkabilov]

If you are considering few other cases, it would be necessary to consider which parameters to be optimized and which one to be left out.

If params_to_optimize = [1>0 1<0 1>0]; then k1 and k_3 are to be considered and k_2 is to be kept const.

If params_to_optimize = [1>0 1>0 1<0]; then k1 and k_2 are to be considered and k_3 is to be kept const.

Etc. .. etc.

That is the intitial and quick thought here.

[Rebeca Hannah Oliveira]

Thanks! The final application is for larger modela (>30 parameters being estimated), so I was trying to find a work-around for this. Please let me know if you have any other ideas. Thanks!