Clear Filters
Clear Filters

How to fix issue with error using surf: Z must be a matrix, not a scalar or vector

4 views (last 30 days)
How can I fix the code to avoid the error: error using surf: Z must be a matrix, not a scalar or vector
I tried converting it to matrix format using cell2mat, but then I can't refer to the correct position of terms to plot on the surface. Is there an alternative way?
numODEs =3;
numParams=3;
% Parameters
parameters = rand(1, numParams); % Initialize parameters with random values
% Initial conditions
initial_conditions = rand(1, numODEs);
% Time span for simulation
tspan = [0 10];
% Set variation range based on a percentage of original parameters
min_variation_factor = 0.1;
max_variation_factor = 10;
variation_range = linspace(min_variation_factor, max_variation_factor, 5);
% Store results
balance_results = cell(numParams, numParams);
% Explore parameter balance for all combinations
for i = 1:numParams
for j = 1:numParams
% Perturb parameters based on variation_range
perturbed_params = parameters;
perturbed_params(i) = perturbed_params(i) * variation_range(1);
perturbed_params(j) = perturbed_params(j) * variation_range(2);
% Solve ODE system with perturbed parameters using ode15s
[~, Y_perturbed] = ode15s(@ode_system, tspan, initial_conditions, [], perturbed_params);
% Check if the solution compensates for variations
balance_results{i, j} = max(Y_perturbed, [], 1); % Adjust based on your specific measure
end
end
balance_results{1, 3}(:, 1)
ans = 0.8496
% Visualize parameter balance
figure;
for k = 1:numODEs
subplot(ceil(sqrt(numODEs)), ceil(sqrt(numODEs)), k);
colormap(parula(numParams^2));
for i = 1:numParams
for j = 1:numParams
surf(variation_range(1), variation_range(2), balance_results{i, j}(:, k)', 'FaceColor', 'interp', 'EdgeColor', 'none');
hold on;
end
end
xlabel('Variation in Parameter 1');
ylabel('Variation in Parameter 2');
zlabel(['Max Output for ODE ', num2str(k)]);
title(['Balance Exploration for ODE ', num2str(k)]);
view(-45, 30); % Adjust the view for better visibility
hold off;
end
Error using surf
Z must be a matrix, not a scalar or vector.
function dydt = ode_system(t, y, parameters)
% ODE system representing binding reactions (modify based on your specific system)
dydt = zeros(size(y));
% Parameters
k1 = parameters(1);
k2 = parameters(2);
k3 = parameters(3);
dydt(1) = -k1 * y(1) * y(2); % Example binding reaction
dydt(2) = k1 * y(1) * y(2) - k2 * y(2); % Example binding and unbinding reaction
dydt(3) = k2 * y(2) - k3 * y(3); % Example unbinding reaction
% Define the ODEs (modify based on your specific system)
% Example: dydt(1) = -parameters(1) * y(1) * y(2);
% Modify the code based on the structure of your ODE system
end
  2 Comments
Dyuman Joshi
Dyuman Joshi on 20 Nov 2023
How do you plan to plot a surface where the data you have is scalar? You can't even plot a line with a single point.
What are you trying to do? What is the objective here?
Rebeca Hannah Oliveira
Rebeca Hannah Oliveira on 20 Nov 2023
Edited: Rebeca Hannah Oliveira on 20 Nov 2023
Thank you! My goal is to get something like in the code below (which was specific for three ODEs and 3 parameters), but generalized to any number n of ODEs and parameters p. I want to generate n figures, each with p subplots that are surface plots showing what happens to the output of the ODEs when two parameters are varied within a range.
% function explore_parameter_balance_binding_reactions()
% Number of ODEs and parameters
numODEs = 3;
numParams = 3;
% Parameters
parameters = rand(1, numParams); % Initialize parameters with random values
% Initial conditions
initial_conditions = rand(1, numODEs);
% Time span for simulation
tspan = [0 10];
% Set variation range based on a percentage of original parameters
min_variation_factor = 0.1;
max_variation_factor = 10;
variation_range = linspace(min_variation_factor, max_variation_factor, 5);
% Store results
balance_results_1_2 = zeros(length(variation_range), length(variation_range), numODEs);
balance_results_1_3 = zeros(length(variation_range), length(variation_range), numODEs);
balance_results_2_3 = zeros(length(variation_range), length(variation_range), numODEs);
% Explore parameter balance for all combinations
for i = 1:length(variation_range)
for j = 1:length(variation_range)
% Perturb parameters based on variation_range
parameters_perturbed_1_2 = parameters .* [variation_range(i), variation_range(j), 1];
parameters_perturbed_1_3 = parameters .* [variation_range(i), 1, variation_range(j)];
parameters_perturbed_2_3 = parameters .* [1, variation_range(i), variation_range(j)];
% Solve ODE system with perturbed parameters using ode15s
[~, Y_perturbed_1_2] = ode15s(@ode_system, tspan, initial_conditions, [], parameters_perturbed_1_2);
[~, Y_perturbed_1_3] = ode15s(@ode_system, tspan, initial_conditions, [], parameters_perturbed_1_3);
[~, Y_perturbed_2_3] = ode15s(@ode_system, tspan, initial_conditions, [], parameters_perturbed_2_3);
% Check if the solution compensates for variations
balance_results_1_2(i, j, :) = max(Y_perturbed_1_2, [], 1); % Adjust based on your specific measure
balance_results_1_3(i, j, :) = max(Y_perturbed_1_3, [], 1);
balance_results_2_3(i, j, :) = max(Y_perturbed_2_3, [], 1);
end
end
% Visualize parameter balance for parameters 1 and 2
figure;
for k = 1:numODEs
subplot(ceil(sqrt(numODEs)), ceil(sqrt(numODEs)), k);
surf(variation_range, variation_range, squeeze(balance_results_1_2(:, :, k))');
xlabel('Variation in Parameter 1');
ylabel('Variation in Parameter 2');
zlabel(['Max Output for ODE ', num2str(k)]);
title(['Balance Exploration for ODE ', num2str(k), ' (Params 1 and 2)']);
end
sgtitle('Parameter Balance Exploration for Binding Reactions (Params 1 and 2)');
% Visualize parameter balance for parameters 1 and 3
figure;
for k = 1:numODEs
subplot(ceil(sqrt(numODEs)), ceil(sqrt(numODEs)), k);
surf(variation_range, variation_range, squeeze(balance_results_1_3(:, :, k))');
xlabel('Variation in Parameter 1');
ylabel('Variation in Parameter 3');
zlabel(['Max Output for ODE ', num2str(k)]);
title(['Balance Exploration for ODE ', num2str(k), ' (Params 1 and 3)']);
end
sgtitle('Parameter Balance Exploration for Binding Reactions (Params 1 and 3)');
% Visualize parameter balance for parameters 2 and 3
figure;
for k = 1:numODEs
subplot(ceil(sqrt(numODEs)), ceil(sqrt(numODEs)), k);
surf(variation_range, variation_range, squeeze(balance_results_2_3(:, :, k))');
xlabel('Variation in Parameter 2');
ylabel('Variation in Parameter 3');
zlabel(['Max Output for ODE ', num2str(k)]);
title(['Balance Exploration for ODE ', num2str(k), ' (Params 2 and 3)']);
end
sgtitle('Parameter Balance Exploration for Binding Reactions (Params 2 and 3)');
% end
function dydt = ode_system(t, y, parameters)
% ODE system representing binding reactions (modify based on your specific system)
dydt = zeros(size(y));
% Parameters
k1 = parameters(1);
k2 = parameters(2);
k3 = parameters(3);
% Define the ODEs (modify based on your specific system)
dydt(1) = -k1 * y(1) * y(2); % Example binding reaction
dydt(2) = k1 * y(1) * y(2) - k2 * y(2); % Example binding and unbinding reaction
dydt(3) = k2 * y(2) - k3 * y(3); % Example unbinding reaction
end

Sign in to comment.

Answers (1)

Pratyush
Pratyush on 20 Nov 2023
Hi Rebecca,
I understand that you are getting an error using the "surf" function in your script.
You could use "meshgrid" function to resolve the error. To fix this, you need to create a meshgrid of "variation_range(1)" and "variation_range(2)" to use as the X and Y inputs for the "surf" function. Then you can use the "meshgrid" function to create the grid and then plot the surfaces accordingly. Here's how you can modify your code:
% Visualize parameter balance
figure;
for k = 1:numODEs
subplot(ceil(sqrt(numODEs)), ceil(sqrt(numODEs)), k);
colormap(parula(numParams^2));
[X, Y] = meshgrid(variation_range, variation_range);
for i = 1:numParams
for j = 1:numParams
Z = balance_results{i, j}(:, k)';
surf(X, Y, reshape(Z, size(X)), 'FaceColor', 'interp', 'EdgeColor', 'none');
hold on;
end
end
xlabel('Variation in Parameter 1');
ylabel('Variation in Parameter 2');
zlabel(['Max Output for ODE ', num2str(k)]);
title(['Balance Exploration for ODE ', num2str(k)]);
view(-45, 30); % Adjust the view for better visibility
hold off;
end
I hope this should resolve the error.
  2 Comments
Dyuman Joshi
Dyuman Joshi on 20 Nov 2023
There's an error when I try to run the code with your suggestion -
numODEs =3;
numParams=3;
% Parameters
parameters = rand(1, numParams); % Initialize parameters with random values
% Initial conditions
initial_conditions = rand(1, numODEs);
% Time span for simulation
tspan = [0 10];
% Set variation range based on a percentage of original parameters
min_variation_factor = 0.1;
max_variation_factor = 10;
variation_range = linspace(min_variation_factor, max_variation_factor, 5);
% Store results
balance_results = cell(numParams, numParams);
% Explore parameter balance for all combinations
for i = 1:numParams
for j = 1:numParams
% Perturb parameters based on variation_range
perturbed_params = parameters;
perturbed_params(i) = perturbed_params(i) * variation_range(1);
perturbed_params(j) = perturbed_params(j) * variation_range(2);
% Solve ODE system with perturbed parameters using ode15s
[~, Y_perturbed] = ode15s(@ode_system, tspan, initial_conditions, [], perturbed_params);
% Check if the solution compensates for variations
balance_results{i, j} = max(Y_perturbed, [], 1); % Adjust based on your specific measure
end
end
balance_results{1, 3}(:, 1)
ans = 0.3748
% Visualize parameter balance
figure;
for k = 1:numODEs
subplot(ceil(sqrt(numODEs)), ceil(sqrt(numODEs)), k);
colormap(parula(numParams^2));
[X, Y] = meshgrid(variation_range, variation_range);
for i = 1:numParams
for j = 1:numParams
Z = balance_results{i, j}(:, k)';
surf(X, Y, reshape(Z, size(X)), 'FaceColor', 'interp', 'EdgeColor', 'none');
hold on;
end
end
xlabel('Variation in Parameter 1');
ylabel('Variation in Parameter 2');
zlabel(['Max Output for ODE ', num2str(k)]);
title(['Balance Exploration for ODE ', num2str(k)]);
view(-45, 30); % Adjust the view for better visibility
hold off;
end
Error using reshape
Number of elements must not change. Use [] as one of the size inputs to automatically calculate the appropriate size for that dimension.
function dydt = ode_system(t, y, parameters)
% ODE system representing binding reactions (modify based on your specific system)
dydt = zeros(size(y));
% Parameters
k1 = parameters(1);
k2 = parameters(2);
k3 = parameters(3);
dydt(1) = -k1 * y(1) * y(2); % Example binding reaction
dydt(2) = k1 * y(1) * y(2) - k2 * y(2); % Example binding and unbinding reaction
dydt(3) = k2 * y(2) - k3 * y(3); % Example unbinding reaction
% Define the ODEs (modify based on your specific system)
% Example: dydt(1) = -parameters(1) * y(1) * y(2);
% Modify the code based on the structure of your ODE system
end

Sign in to comment.

Categories

Find more on Programming in Help Center and File Exchange

Products


Release

R2023b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!