getting error in code when running

CheckKFmax_1()

R*peak: 65.4697
T-peak: 4.5
50% form peak to ss : 57.7799
R* Steady state: 50.0901
T-50 : 53
R*peak/R*ss: 1.307
Δ (R*peak-R*ss) : 15.3795
R*peak: 63.795
T-peak: 3.7
50% form peak to ss : 48.6249
R* Steady state: 33.4548
T-50 : 39.1
R*peak/R*ss: 1.9069
Δ (R*peak-R*ss) : 30.3402
R*peak: 62.4078
T-peak: 3.2
50% form peak to ss : 43.7608
R* Steady state: 25.1139
T-50 : 31.3
R*peak/R*ss: 2.485
Δ (R*peak-R*ss) : 37.2938
R*peak: 61.2429
T-peak: 2.9
50% form peak to ss : 40.6736
R* Steady state: 20.1044
T-50 : 26.2
R*peak/R*ss: 3.0462
Δ (R*peak-R*ss) : 41.1385
R*peak: 60.17
T-peak: 2.7
50% form peak to ss : 38.4641
R* Steady state: 16.7582
T-50 : 22.6
R*peak/R*ss: 3.5905
Δ (R*peak-R*ss) : 43.4118
R*peak: 79.1276
T-peak: 4.7
50% form peak to ss : 72.9374
R* Steady state: 66.7472
T-50 : 56.7
R*peak/R*ss: 1.1855
Δ (R*peak-R*ss) : 12.3804
R*peak: 77.8834
T-peak: 3.7
50% form peak to ss : 64.0101
R* Steady state: 50.1368
T-50 : 46.4
R*peak/R*ss: 1.5534
Δ (R*peak-R*ss) : 27.7466
R*peak: 76.843
T-peak: 3.2
50% form peak to ss : 58.5003
R* Steady state: 40.1576
T-50 : 39.4
R*peak/R*ss: 1.9135
Δ (R*peak-R*ss) : 36.6854
R*peak: 75.9321
T-peak: 3
50% form peak to ss : 54.7034
R* Steady state: 33.4747
T-50 : 34
R*peak/R*ss: 2.2683
Δ (R*peak-R*ss) : 42.4574
R*peak: 75.0963
T-peak: 2.8
50% form peak to ss : 51.9009
R* Steady state: 28.7056
T-50 : 30.1
R*peak/R*ss: 2.6161
Δ (R*peak-R*ss) : 46.3908
R*peak: 85.0465
T-peak: 4.7
50% form peak to ss : 80.0572
R* Steady state: 75.0679
T-50 : 58.7
R*peak/R*ss: 1.1329
Δ (R*peak-R*ss) : 9.9786
R*peak: 84.0739
T-peak: 3.6
50% form peak to ss : 72.099
R* Steady state: 60.1241
T-50 : 50.8
R*peak/R*ss: 1.3983
Δ (R*peak-R*ss) : 23.9498
R*peak: 83.2706
T-peak: 3.2
50% form peak to ss : 66.7106
R* Steady state: 50.1506
T-50 : 44.5
R*peak/R*ss: 1.6604
Δ (R*peak-R*ss) : 33.12
R*peak: 82.5506
T-peak: 2.9
50% form peak to ss : 62.7819
R* Steady state: 43.0133
T-50 : 39.7
R*peak/R*ss: 1.9192
Δ (R*peak-R*ss) : 39.5374
R*peak: 81.889
T-peak: 2.7
50% form peak to ss : 59.7715
R* Steady state: 37.6539
T-50 : 35.7
R*peak/R*ss: 2.1748
Δ (R*peak-R*ss) : 44.2351
R*peak: 88.3479
T-peak: 4.7
50% form peak to ss : 84.2016
R* Steady state: 80.0553
T-50 : 60.2
R*peak/R*ss: 1.1036
Δ (R*peak-R*ss) : 8.2926
R*peak: 87.5604
T-peak: 3.8
50% form peak to ss : 77.1745
R* Steady state: 66.7885
T-50 : 53.1
R*peak/R*ss: 1.311
Δ (R*peak-R*ss) : 20.772
R*peak: 86.9035
T-peak: 3.2
50% form peak to ss : 72.0973
R* Steady state: 57.291
T-50 : 47.9
R*peak/R*ss: 1.5169
Δ (R*peak-R*ss) : 29.6124
R*peak: 86.3115
T-peak: 3
50% form peak to ss : 68.2354
R* Steady state: 50.1594
T-50 : 43.3
R*peak/R*ss: 1.7207
Δ (R*peak-R*ss) : 36.1521
R*peak: 85.7664
T-peak: 2.7
50% form peak to ss : 65.1864
R* Steady state: 44.6064
T-50 : 39.8
R*peak/R*ss: 1.9227
Δ (R*peak-R*ss) : 41.16
R*peak: 90.4558
T-peak: 4.6
50% form peak to ss : 86.9198
R* Steady state: 83.3838
T-50 : 61.4
R*peak/R*ss: 1.0848
Δ (R*peak-R*ss) : 7.072
R*peak: 89.7938
T-peak: 3.8
50% form peak to ss : 80.6666
R* Steady state: 71.5394
T-50 : 55
R*peak/R*ss: 1.2552
Δ (R*peak-R*ss) : 18.2545
R*peak: 89.2374
T-peak: 3.3
50% form peak to ss : 75.9399
R* Steady state: 62.6424
T-50 : 50.3
R*peak/R*ss: 1.4246
Δ (R*peak-R*ss) : 26.595
R*peak: 88.7405
T-peak: 2.9
50% form peak to ss : 72.2258
R* Steady state: 55.711
T-50 : 46.4
R*peak/R*ss: 1.5929
Δ (R*peak-R*ss) : 33.0294
R*peak: 88.2778
T-peak: 2.7
50% form peak to ss : 69.2131
R* Steady state: 50.1483
T-50 : 42.9
R*peak/R*ss: 1.7603
Δ (R*peak-R*ss) : 38.1295

Error using .
Unrecognized table variable name 'Kf1Max'.

Error in solution>CheckKFmax_1 (line 86)
Kf1Max_values = unique(combined_table.Kf1Max);

function CheckKFmax_1()
    % Parameters
    timespan = 0:0.1:500;  % Time vector
    Receptor_concentration = 100;  % Concentration of Receptor
    
    % Initial conditions
    C_LigandReceptor_0 = 0;  % Initial concentration of the ligand-receptor complex
    % Define different ranges of Kf1Max values
    Kf1Max_values = [20, 40, 60, 80, 100];
    Kb1Max_values = [20, 40, 60, 80, 100];
    % Create a cell array to store tables for different Kf1Max and Kb1Max ranges
    tables_cell = cell(length(Kf1Max_values), length(Kb1Max_values));
    for i = 1:length(Kf1Max_values)
        Kf1Max = Kf1Max_values(i);
        for j = 1:length(Kb1Max_values)
            Kb1Max = Kb1Max_values(j);
            % Time-dependent forward reaction rate constants
            t_start = 5;   % Start time of receptor pulse
            t_end = 500;    % Duration of receptor pulse
            Kf1Min = 1;
            Kf1Tau_on = -1;
            Kf1Tau_off = -1;
  
            % Time-dependent backward reaction rate constants
            Kb1Min = 10;
            Kb1Tau_on =  -0.01;
            Kb1Tau_off = -0.01;
            % Define anonymous functions for forward and backward reaction rates
            kf_1 = @(t) calculate_kf(t, t_start, t_end, Kf1Max, Kf1Min, Kf1Tau_on, Kf1Tau_off);
            kb_1 = @(t) calculate_kb(t, t_start, t_end, Kb1Max, Kb1Min, Kb1Tau_on, Kb1Tau_off);
            % Solve the ODE system
            [t, y] = ode45(@(t, y) ode_LR(t, y, kf_1, kb_1), timespan, [Receptor_concentration; C_LigandReceptor_0]);
            % Extract the concentrations
            receptor_concentration = y(:, 1);
            ligand_Receptor_concentration = y(:, 2);
            % Displaying the peak value and its corresponding time
            [peak_value, peak_time_idx] = max(ligand_Receptor_concentration);
            time_to_peak = t(peak_time_idx);
            % Calculate the time between start and peak time
            time_to_peak = time_to_peak - t_start;
            % Find steady state and time to approach 50% of the steady state
            steady_state = mean(ligand_Receptor_concentration(end-50:end)); % Considering the last few points for steady state
            % Calculate T-50: Time to reach 50% of the peak value
            half_peak_value =  (steady_state + peak_value)/2;
            % Find the indices where the ligand receptor concentration is closest to half the peak value
            [~, idx_50_percent] = min(abs(ligand_Receptor_concentration(peak_time_idx:end) - half_peak_value));
            % Get the time corresponding to the closest index
            time_to_50_percent = t(idx_50_percent)- time_to_peak;
            % Find the ratio of the R*peak/R* steady state 
            peak_to_steady_state_ratio = peak_value / steady_state;
            % Find the ratio of the R*peak/R* steady state 
            Delta = peak_value - steady_state;
            disp(['R*peak: ', num2str(peak_value)]);
            disp(['T-peak: ', num2str(time_to_peak)]);
            disp(['50% form peak to ss : ', num2str(half_peak_value)]);
            disp(['R* Steady state: ', num2str(steady_state)]);
            disp(['T-50 : ', num2str(time_to_50_percent)]);
            disp(['R*peak/R*ss: ', num2str(peak_to_steady_state_ratio)]);
            disp(['Δ (R*peak-R*ss) : ', num2str(Delta)]);
            R_peak_values(i) = peak_value;
            T_peak_values(i) = time_to_peak;
            R_steady_state_values(i) = steady_state;
            T_50_values(i) = time_to_50_percent;
            R_peak_R_ss_values(i) = peak_to_steady_state_ratio;
            Delta_values(i) = Delta;
            % Create a table to combine all results
            data_table = table(R_peak_values', T_peak_values', R_steady_state_values', T_50_values', R_peak_R_ss_values', Delta_values', ...
                'VariableNames', {'R_peak', 'T_peak', 'R_steady_state', 'T_50', 'R_peak_R_ss', 'Delta'});
            % Store the table in the cell array for each Kf1Max range
            tables_cell{i, j} = data_table;
            % Generate a unique filename for each combination of Kf1Max and Kb1Max values
            csv_filename = sprintf('interaction_data_Kf1Max_%d_Kb1Max_%d.csv', Kf1Max, Kb1Max);
            writetable(data_table, csv_filename);
        end
    end
    % Combine all results into a single table
    combined_table = vertcat(tables_cell{:});
    % Generate a unique filename for the combined table
    combined_csv_filename = 'combined_interaction_data_all.csv';
    writetable(combined_table, combined_csv_filename);
    % Load the combined interaction data CSV file
    combined_csv_filename = 'combined_interaction_data_all.csv';
    combined_table = readtable(combined_csv_filename);
    % Extract Kf1Max and Kb1Max values from the combined table
    Kf1Max_values = unique(combined_table.Kf1Max);
    Kb1Max_values = unique(combined_table.Kb1Max);
    % Initialize matrices to store Kf1Max, Kb1Max, and R*peak values
    Kf1Max_matrix = zeros(length(Kb1Max_values), length(Kf1Max_values));
    Kb1Max_matrix = zeros(length(Kb1Max_values), length(Kf1Max_values));
    R_peak_matrix = zeros(length(Kb1Max_values), length(Kf1Max_values));
    % Populate matrices with values from the combined table
    for i = 1:length(Kf1Max_values)
        for j = 1:length(Kb1Max_values)
            idx = find(combined_table.Kf1Max == Kf1Max_values(i) & combined_table.Kb1Max == Kb1Max_values(j));
            if ~isempty(idx)
                Kf1Max_matrix(j, i) = Kf1Max_values(i);
                Kb1Max_matrix(j, i) = Kb1Max_values(j);
                R_peak_matrix(j, i) = combined_table.R_peak(idx);
            end
        end
    end
    % Create a 3D fsurf plot
    figure;
    fsurf(Kf1Max_matrix, Kb1Max_matrix, R_peak_matrix);
    xlabel('Kf1Max');
    ylabel('Kb1Max');
    zlabel('R*peak');
    title('3D fsurf Plot for Combined R*peak Values');
end
function kf_1 = calculate_kf(t, t_start, t_end, Kf1Max, Kf1Min, Kf1Tau_on, Kf1Tau_off)
    % Initialize kf_1
    kf_1 = zeros(size(t));
    
    % Initial phase: Constant value Kf1Min
    initial_phase = t < t_start;
    kf_1(initial_phase) = Kf1Min;
    
    % Increasing phase: Exponential growth from Kf1Min to Kf1Max
    increasing_phase = t >= t_start & t < t_end;
    kf_1(increasing_phase) = Kf1Max - (Kf1Max - Kf1Min) * exp(Kf1Tau_on * (t(increasing_phase) - t_start));
    
    % Decreasing phase: Exponential decay from Kf1Max to Kf1Min
    kf_end = Kf1Max - (Kf1Max - Kf1Min) * exp(Kf1Tau_on * (t_end - t_start));
    decreasing_phase = t >= t_end;
    kf_1(decreasing_phase) = Kf1Min + (kf_end - Kf1Min) * exp(Kf1Tau_off * (t(decreasing_phase) - t_end));
end
function kb_1 = calculate_kb(t, t_start, t_end, Kb1Max, Kb1Min, Kb1Tau_on, Kb1Tau_off)
    % Initialize kb_1
    kb_1 = zeros(size(t));
    
    % Initial phase: Constant value Kb1Min
    initial_phase = t < t_start;
    kb_1(initial_phase) = Kb1Min;
    
    % Increasing phase: Exponential growth from Kb1Min to Kb1Max
    increasing_phase = t >= t_start & t < t_end;
    kb_1(increasing_phase) = Kb1Max - (Kb1Max - Kb1Min) * exp(Kb1Tau_on * (t(increasing_phase) - t_start));
    
    % Decreasing phase: Exponential decay from Kb1Max to Kb1Min
    kb_end = Kb1Max - (Kb1Max - Kb1Min) * exp(Kb1Tau_on * (t_end - t_start)); 
    decreasing_phase = t >= t_end;
    kb_1(decreasing_phase) = Kb1Min + (kb_end - Kb1Min) * exp(Kb1Tau_off * (t(decreasing_phase) - t_end));
end
function dydt = ode_LR(t, y, kf_1, kb_1)
    % Unpack the concentrations
    Receptor_concentration = y(1);
    C_LigandReceptor = y(2);
    
    % Define the ODE system
    dReceptor_dt = -kf_1(t) * Receptor_concentration + kb_1(t) * C_LigandReceptor;
    d_CLigandReceptor_dt = kf_1(t) * Receptor_concentration - kb_1(t) * C_LigandReceptor;
    % Pack the derivatives into a column vector
    dydt = [dReceptor_dt; d_CLigandReceptor_dt];
end

Answer 1

[Harald]

Hi,

after running your code, combined_table does not contain anything called Kf1Max or even remotely similar to it.

Best wishes,

Harald

[Ehtisham]

i want to get a 3D mesh plot of R*peak values for all the iteration in loop.

[Image Analyst]

OK, fine - the values are already there in the table. It looks like Kf1Max is not even involved. Do you have an (x,y) coordinate for each value of R_peak? If so, then use surf

[Harald]

@Ehtisham, it is of course always good to know where you ultimately want to go. However, I don't think this resolves the imminent problem: you are trying to extract a variable from a table that isn't there.

So, as Image Analyst suggests, you will need to obtain that variable from "somewhere" else.