CDF of log pearson
5 ビュー (過去 30 日間)
古いコメントを表示
Can somebody explain why when I use this code, I get CDF as negative and decreasing function (magnitude)? It should be non decreasing function
where Q include the data given below
23.81
33.98
62.01
140.45
184.91
257.97
325.64
410.59
543.68
I have integrated the pdf to calculate cdf of log pearson type III distribution
MATLAB CODE
Q = load('F:\DISCHARGE data.txt')
n = length(Q);
Z = log(Q);
mu = mean(Z);%%mean
Sigma = std(Z);%%Standard deviattion%
Skew = sum(n*((Z-mu).^3))/((Sigma.^3)*(n-1)*(n-2));%%Skewness%%
kurt = (((n*(n+1))/((n-1)*(n-2)*(n-3)))*(sum(((Z-mu)/Sigma).^4)))-((3*((n-1).^2))/((n-2)*(n-3)));
alpha_p = 4/(Skew^2);%%Shape_alpha_b%%
beta_b = (Sigma)/(alpha_p^0.5); %%scale_beta_1bya%%
y_a = mu-(Sigma*(alpha_p^0.5));%%location_y_m%%
syms e
fun = (1/(e*gamma(alpha_p).*abs(beta_b))).*(((log(e)-y_a)/beta_b).^(alpha_p - 1)).*(exp(-((log(e)-y_a)/beta_b)));
wer = vpa(simplify(int(fun,e)))
CDF = round(double(subs(wer,Q)),4)
OUTPUT
-0.9674
-0.9175
-0.7612
-0.4678
-0.3738
-0.2743
-0.2158
-0.1669
-0.1196
0 件のコメント
採用された回答
Paul
2021 年 7 月 10 日
編集済み: Paul
2021 年 7 月 11 日
Not familiar with the disribution, but can offer the following about this code:
% data
Q=[23.81
33.98
62.01
140.45
184.91
257.97
325.64
410.59
543.68];
% parameters
n = length(Q);
Z = log(Q);
mu = mean(Z);%%mean
Sigma = std(Z);%%Standard deviattion%
Skew = sum(n*((Z-mu).^3))/((Sigma.^3)*(n-1)*(n-2));%%Skewness%%
kurt = (((n*(n+1))/((n-1)*(n-2)*(n-3)))*(sum(((Z-mu)/Sigma).^4)))-((3*((n-1).^2))/((n-2)*(n-3)));
alpha_p = 4/(Skew^2);%%Shape_alpha_b%%
beta_b = (Sigma)/(alpha_p^0.5); %%scale_beta_1bya%%
y_a = mu-(Sigma*(alpha_p^0.5));%%location_y_m%%
% pdf
syms e positive
fun(e) = (1/(e*gamma(alpha_p).*abs(beta_b))).*(((log(e)-y_a)/beta_b).^(alpha_p - 1)).*(exp(-((log(e)-y_a)/beta_b)));
The support of the pdf is not given. Note that for e less than 3 fun(e) is complex, and f(3) ~ 0
double(fun(2))
double(fun(3))
So we'll assume the support of fun is 3 < e < inf.
The CDF is the integral of the pdf. Note that we must integrate from 3 to q to get P(3 < Q < q). Compare this to the original code
syms q positive
F(q) = simplify(int(fun(e),e,3,q));
Plot the CDF:
fplot(F(q),[3 1000])
That looks like a CDF (note that F(q) = 0 for q < 3). No idea if it's the CDF you're expecting.
Looks like it's very close to, if not exactly, the incomplete gamma function with support x >= exp(y_a)
hold on;
x = 3:10:1000;
plot(x,gammainc((log(x)-y_a)/beta_b,alpha_p),'ro')
その他の回答 (0 件)
参考
カテゴリ
Help Center および File Exchange で Probability Distributions についてさらに検索
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!