I want to create a QQ & PP plot using L-Moments to fit the GEV but get this error: unable to perform assignment b/c the left & right sides have a different number of elements

Here is my code:

% Step A: Read data and extract columns
data = readmatrix('Wabash_River_MaxFlow.xls');

Error using readmatrix
Unable to find or open 'Wabash_River_MaxFlow.xls'. Check the path and filename or file permissions.

wateryear = data(:, 1);
streamflow = data(:, 2);

% Step B: Compute L-Moments
n = length(data);
b0 = mean(data);
run_sum = 0;
for i = 1:(n-1)
    run_sum = run_sum + (i+1-1)*data(i+1);
end
b1 = run_sum/(n*(n-1));

run_sum = 0;
for i = 2:(n-1)
    run_sum = run_sum + (i+1-1)*(i+1-2)*data(i+1);
end
b2 = run_sum/(n*(n-1)*(n-2));

run_sum = 0;
for i = 3:(n-1)
    run_sum = run_sum + (i+1-1)*(i+1-2)*(i+1-3)*data(i+1);
end
b3 = run_sum/(n*(n-1)*(n-2)*(n-3));

run_sum = 0;
for i = 4:(n-1)
    run_sum = run_sum + (i+1-1)*(i+1-2)*(i+1-3)*(i+1-4)*data(i+1);
end
b4 = run_sum/(n*(n-1)*(n-2)*(n-3)*(n-4));

l1 = b0;
l2 = 2*b1-b0;
l3 = 6*b2-6*b1+b0;
l4 = 20*b3-30*b2+12*b1-b0;

L_CV = l2/l1;
L_SK = l3/l2;
L_KU = l4/l2;

% Step C: Generate QQ plot
X_plot = zeros(1, 80);
pp_plot = linspace(1/81, 1, 80);  % Update: generate pp_plot directly

for i = 1:80
    X_plot(i) = psi + alpha*(1-(-1*log(pp_plot(i)))^k)/k;
end

figure;
scatter(pp_plot, streamflow);
hold on;
plot(pp_plot, X_plot, '-r');
xlabel('Non-Exceedance Prob.');
ylabel('Observations');
title('QQ Plot');
legend('Observations', 'GEV Fit');

% Step D: Generate PP plot
pp_data = linspace(1/(n+1), n/(n+1), n);  % Update: generate pp_data directly
pp_theory = psi + alpha*(1-(-1*log(pp_data))^k)/k;

figure;
plot(linspace(0, 1, 100), linspace(0, 1, 100), '-r');
hold on;
scatter(pp_theory, pp_data);
xlabel('Theoretical Non-Exceedance Prob.');
ylabel('Fitted Non-Exceedance Prob.');
title('PP Plot');
legend('1:1 Line', 'GEV Fit');

Here is the error I am getting:

Unable to perform assignment because the left and right sides have a different
number of elements.

Error in CEE204_HW2_Problem1 (line 47)
    X_plot(i) = psi + alpha*(1-(-1*log(pp_plot(i)))^k)/k;

Here is what my data looks like:

Answer 1

[Jeff Miller]

Assign the problematic result to a new variable so that you can see what elements it has, like this:

 temp_var = psi + alpha*(1-(-1*log(pp_plot(i)))^k)/k;

Probably you will see that temp_var does not have the number of elements you expected. For example, maybe instead you need:

 X_plot(i) = psi + alpha*(1-(-1*log(pp_plot(i))).^k)/k;