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
4 Ansichten (letzte 30 Tage)
Ältere Kommentare anzeigen
Here is my code:
% Step A: Read data and extract columns
data = readmatrix('Wabash_River_MaxFlow.xls');
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:
0 Kommentare
Antworten (1)
Jeff Miller
am 6 Apr. 2023
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;
0 Kommentare
Siehe auch
Kategorien
Mehr zu Generalized Extreme Value Distribution finden Sie in Help Center und File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!