I have fitted a distribution curve using gamma probability density distribution with 3 parameters. The following are my codes:
proc nlin data=worker_distribution
method=marquardt ;
parameters alpha=3 beta=2 lambda=0.01 ;
model percentcount=pdf("gamma",ratio,alpha,beta)*lambda;
output out=fitworker predicted=predworker parms= alpha beta lambda;
run;
How do I create the cdf in SAS EG based on the fitted gamma pdf distribution? Thanks.
> I did use lambda to get the final fitted gamma pdf. Lambda is a constant and we integrate pdf to get cdf , so if I simply multiply lambda by cdf("gamma", x, alpha, beta), I get very small number and it is no way near 1. What am I missing here?
You might be missing that there are two quantities being estimated. You use lambda to estimate the model for the percentCount variable. The predicted mean of percentCount at the ratio r is
predworker = lambda*(pdf("gamma", r, alpha, beta)
Notice that predworker is not a density. It does not integrate to unity. It will integrate to lambda, which for your data is approximately 0.01, which is "a very small number."
The following program simulates some data, then uses PROC NLIN to estimate the parameters. I then graph the PDF and the predicted mean for PredWorker. The PDF integrates to unity. The predicted model does not integrate to unity.
data Have;
call streaminit(1);
alpha=3;
beta=0.2;
lambda=0.01;
do i = 1 to 100;
ratio = rand("uniform", 0, 3);
percentcount=pdf("gamma",ratio,alpha,beta)*lambda + rand("Normal",0, 0.001);
output;
end;
drop alpha beta lambda;
run;
title "Observed PercentCount per Ratio";
proc sgplot data=Have;
scatter x=ratio y=percentcount;
run;
proc nlin data=have method=marquardt ;
parameters alpha=3 beta=2 lambda=0.01 ;
model percentcount=pdf("gamma",ratio,alpha,beta)*lambda;
output out=fitworker predicted=predworker parms= alpha beta lambda;
run;
proc sort data=fitworker; by ratio; run;
proc sgplot data=fitworker;
scatter x=ratio y=percentcount;
series x=ratio y=predworker;
run;
data Gamma;
set fitworker(obs=1 keep=alpha beta lambda);
do x = 0.05 to 3 by 0.05;
pdf = PDF("gamma", x, alpha, beta);
cdf = CDF("gamma", x, alpha, beta);
pred = pdf*lambda;
output;
end;
run;
title "PDF('Gamma') and Predicted Values";
proc sgplot data=Gamma;
series x=x y=PDF;
series x=x y=pred;
yaxis grid;
xaxis grid;
run;
If this does not answer your question, you'll have to post some sample data and the code you are using to integrate whatever function you are integrating.
Calling @Rick_SAS
If you have the estimates for the alpha and beta parameters in the gamma distribution, you can use the PDF function to visualize the density and use the CDF function to visualize the cumulative probability, as discussed in the article, "Four essential functions for statistical programmers."
For example, suppose the estimate of the shape parameter is 3.2 and the estimate for the scale parameter is 2.1. Then the following DATA step evaluates the PD and CDF at a set of values in [0, 20] and plots the result:
data Gamma;
alpha = 3.2; /* shape parameter */
beta = 2.1; /* scale parameter */
do x = 0.1 to 20 by 0.1;
pdf = PDF("gamma", x, alpha, beta);
cdf = CDF("gamma", x, alpha, beta);
output;
end;
run;
title "PDF and CDF of the Gamma Distribution";
title2 "Shape=3.2; Scale=2.1";
proc sgplot data=Gamma;
series x=x y=PDF;
series x=x y=CDF;
yaxis grid;
xaxis grid;
run;
What are the estimates of alpha, beta, and lambda?
I appreciate your quick responses. Here are the parameters from the fitting procedure and I attached pdf graphs. Once again, thank you for your help.
Alpha=3.1264265313
Beta=0.1878574508
Lambda=0.0096286931
> I did use lambda to get the final fitted gamma pdf. Lambda is a constant and we integrate pdf to get cdf , so if I simply multiply lambda by cdf("gamma", x, alpha, beta), I get very small number and it is no way near 1. What am I missing here?
You might be missing that there are two quantities being estimated. You use lambda to estimate the model for the percentCount variable. The predicted mean of percentCount at the ratio r is
predworker = lambda*(pdf("gamma", r, alpha, beta)
Notice that predworker is not a density. It does not integrate to unity. It will integrate to lambda, which for your data is approximately 0.01, which is "a very small number."
The following program simulates some data, then uses PROC NLIN to estimate the parameters. I then graph the PDF and the predicted mean for PredWorker. The PDF integrates to unity. The predicted model does not integrate to unity.
data Have;
call streaminit(1);
alpha=3;
beta=0.2;
lambda=0.01;
do i = 1 to 100;
ratio = rand("uniform", 0, 3);
percentcount=pdf("gamma",ratio,alpha,beta)*lambda + rand("Normal",0, 0.001);
output;
end;
drop alpha beta lambda;
run;
title "Observed PercentCount per Ratio";
proc sgplot data=Have;
scatter x=ratio y=percentcount;
run;
proc nlin data=have method=marquardt ;
parameters alpha=3 beta=2 lambda=0.01 ;
model percentcount=pdf("gamma",ratio,alpha,beta)*lambda;
output out=fitworker predicted=predworker parms= alpha beta lambda;
run;
proc sort data=fitworker; by ratio; run;
proc sgplot data=fitworker;
scatter x=ratio y=percentcount;
series x=ratio y=predworker;
run;
data Gamma;
set fitworker(obs=1 keep=alpha beta lambda);
do x = 0.05 to 3 by 0.05;
pdf = PDF("gamma", x, alpha, beta);
cdf = CDF("gamma", x, alpha, beta);
pred = pdf*lambda;
output;
end;
run;
title "PDF('Gamma') and Predicted Values";
proc sgplot data=Gamma;
series x=x y=PDF;
series x=x y=pred;
yaxis grid;
xaxis grid;
run;
If this does not answer your question, you'll have to post some sample data and the code you are using to integrate whatever function you are integrating.
Are you ready for the spotlight? We're accepting content ideas for SAS Innovate 2025 to be held May 6-9 in Orlando, FL. The call is open until September 25. Read more here about why you should contribute and what is in it for you!
What’s the difference between SAS Enterprise Guide and SAS Studio? How are they similar? Just ask SAS’ Danny Modlin.
Find more tutorials on the SAS Users YouTube channel.