Inverse Cumulative distribution function

3 visualizaciones (últimos 30 días)
Bashar AlHalaq
Bashar AlHalaq el 30 de Mzo. de 2021
Comentada: Jeff Miller el 6 de Abr. de 2021
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

Respuesta aceptada

Jeff Miller
Jeff Miller el 31 de Mzo. de 2021
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 comentarios
Bashar AlHalaq
Bashar AlHalaq el 6 de Abr. de 2021
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 el 6 de Abr. de 2021
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.

Iniciar sesión para comentar.

Más respuestas (1)

Riley
Riley el 30 de Mzo. de 2021
You can use function ksdensity.

Etiquetas

Community Treasure Hunt

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

Start Hunting!

Translated by