Inverse Cumulative distribution function

6 ビュー (過去 30 日間)
Bashar AlHalaq
Bashar AlHalaq 2021 年 3 月 30 日
コメント済み: Jeff Miller 2021 年 4 月 6 日
Hi,
if i have the following CDF :
F(X)=1/1+b(1+b-(1+b+qx/b))(1+x/b)^-q
How can I find x=inverse(F(x))
with my best reguards

採用された回答

Jeff Miller
Jeff Miller 2021 年 3 月 31 日
Suppose you have a function called thisCDF which computes the CDF for your distribution. Then define
function x = inverseF(p)
inverseFerr = @(x) thisCDF(x) - p;
x = fzero(inverseFerr,0.5); % replace 0.5 with a starting guess for x
end
Now
x = inverseF(0.5) % compute the median, etc
  5 件のコメント
Bashar AlHalaq
Bashar AlHalaq 2021 年 4 月 6 日
I wrote the following code for estimate the parameters b and q:
clc;
clear;
b=0.6;
q=1.5;
T=1;
n=150;
for t=1:T
for i=1:n
x(i)=inverseF(b,q);
xx=sort(x);
pdf=q.*(b.^2+q*xx-xx)./((b.^3+b.^2).*(1+(xx./b)).^(1+q));
F=((1+b-(1+b+q.*xx/b).*(1+xx/b).^-q))/(1+b);
R=1-((1+b-(1+b+(q.*xx./b)).*(1+(xx./b)).^-q))/(1+b);
end
%% 1-MLE
[par1 f]=fsolve(@(S) MLE(xx,S),[1 1]);
b_mle(t)=par1(1); q_mle(t)=par1(2);
f_mle=q_mle.*(b_mle.^2+q_mle*xx-xx)./((b_mle.^3+b_mle.^2).*(1+(xx./b_mle)).^(1+q_mle));
F_mle=((1+b_mle-(1+b_mle+q_mle.*xx/b_mle).*(1+xx/b_mle).^-q_mle))/(1+b_mle);
R_mle=1-(((1+b_mle-(1+b_mle+q_mle.*xx/b_mle).*(1+xx/b_mle).^-q_mle))/(1+b_mle));
%% 2-MOM
[par2 f]=fsolve(@(S) MOM_Method(xx,S),[b q]);
b_mom(t)=par2(1); q_mom(t)=par2(2);
f_mom=q_mom.*(b_mom.^2+q_mom*xx-xx)./((b_mom.^3+b_mom.^2).*(1+(xx./b_mom)).^(1+q_mom));
F_mom=((1+b_mom-(1+b_mom+q_mom.*xx/b_mom).*(1+xx/b_mom).^-q_mom))/(1+b_mom);
R_mom=1-(((1+b_mom-(1+b_mom+q_mom.*xx/b_mom).*(1+xx/b_mom).^-q_mom))/(1+b_mom));
end
%% Plot PDF
figure(1)
plot(pdf,'linewidth',2)
hold on;
plot(f_mle,'linewidth',2)
plot(f_mom,'linewidth',2)
legend('f','f_mle','f_mom')
xlabel('x')
ylabel('f(x)')
title('PDF for DTWPD distrbution')
%% Plot CDF
figure(2)
plot(F,'linewidth',2)
hold on;
plot(F_mle,'linewidth',2)
plot(F_mom,'linewidth',2)
legend('F','F_mle','F_mom')
xlabel('x')
ylabel('F(x)')
title('CDF for DTWPD distrbution')
%% Plot R
figure(3)
plot(R,'linewidth',2)
hold on;
plot(R_mle,'linewidth',2)
plot(R_mom,'linewidth',2)
legend('R','R_mle','R_mom')
xlabel('t')
ylabel('R(t)')
title('Reliability for mixlo distrbution')
with sub function:
1) for generate data:
function x = inverseF(b,q)
u=rand;
inverseFerr = @(x) ((1+b-(1+b+q*x/b)*(1+x/b)^-q))/(1+b);
x = inverseFerr(u);
end
2) mle
function [F]=MLE(xx,S)
n=length(xx);
b=S(1);
q=S(2);
k1=(n*(3*b^2+2*b));
k2=(1+q);
k3=(n/q);
k4=(2*b);
s1=0; s2=0; s3=0; s4=0;
for i=1:n
s1=s1+(1/(b^2+q*xx(i)-xx(i)));
s2=s2+((b^-2)/((1+xx(i))/b));
s3=s3+xx(i)/(b^2+q*xx(i)-xx(i));
s4=s4+log((1+xx(i))/b);
end
F=[k4*s1-k1+k2*s2 k3+s3-s4];
3)mom
function [F]=MOM_Method(xx,S)
m1=mean(xx);
b=S(1);
q=S(2);
n=length(xx);
s1=0;
for i=1:n
s1=s1+xx(i)^2;
end
m2=s1/n;
M1=b^2/(1+b)*(q-1)+2*b/(1+b)*(q-2);
M2=2*b^3/(1+b)*(q-1)*(q-2)+6*b^2/(1+b)*(q-2)*(q-3);
F=[sqrt((M1-m1)^2) sqrt((M2-m2)^2)];
please are these code correct or not لاecause I noticed the curve of the cumulative distribution function stabilizes at the value 0.7 and does not reach (1) . Is the data generation correct or not?
thank you for your attention
Jeff Miller
Jeff Miller 2021 年 4 月 6 日
I don't know what you are trying to compute so I can't say for sure whether your code is correct, but a few things make me suspect it is not correct:
  • The CDF should reach 1 and not stabilize at 0.7. Maybe your CDF function is not correct (e.g., why does it start with 1/1+... which is the same as 1+....?). What distribution are you actually trying to work with?
  • Your inverseF function does not call fzero so it is not getting x values with the cdf of u. Instead, it is computing the CDF of u.

サインインしてコメントする。

その他の回答 (1 件)

Riley
Riley 2021 年 3 月 30 日
You can use function ksdensity.

タグ

Community Treasure Hunt

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

Start Hunting!

Translated by