BookmarkSubscribeRSS Feed
🔒 This topic is solved and locked. Need further help from the community? Please sign in and ask a new question.
DaLack
Calcite | Level 5

I am interested in estimating a distribution given several quantiles (0.05, 0.5, 0.95) and if necessary the mean. The shape should be somehow normal. It exists a useful package called "‘rriskDistributions" in R, explained here

https://cran.r-project.org/web/packages/rriskDistributions/rriskDistributions.pdf

I am looking for an analogue in SAS 9.4. Does anybody have an idea?

Thanks in advance.

1 ACCEPTED SOLUTION

Accepted Solutions
Ksharp
Super User

Check @Rick  's blog 

http://blogs.sas.com/content/iml/2014/06/18/distribution-from-quantiles.html 

 

 

Here is mine 

 

data have;
infile cards truncover expandtabs;
input Quantile Estimate ;
prob=dif(Quantile);
width=dif(Estimate);
cards;
0 80
.01 128
.05 151
.10 162
.25 184
.50 209
.75 236
.90 265
.95 284
.99 333
1 727
;
run;
proc sql;
select prob into : probs separated by ','
from have
where prob is not missing;
select Estimate into : points separated by ','
from have
where Quantile ne 1;
select width into : widths separated by ','
from have
where width is not missing;
select count(*)-1 into : count from have;
quit;
%let n=10000; /* sample size */
data want;
call streaminit(1234);
array points(&count) _temporary_ (&points);
array widths(&count) _temporary_ (&widths);
do i=1 to &n;
n=rand('table',&probs);
x=points{n}+ceil(widths{n}*rand('uniform'));
output;
end;
drop i;
run;
proc univariate data=want;
var x;
histogram x/kernel;
run;

 

 

 

 

Or An fast way is IML code , also written by me . Hope you like it .

 

proc iml;
/*
quantiles of total cholesterol from NHANES study
http://blogs.sas.com/content/iml/2014/06/18/distribution-from-quantiles.html
*/
Quantile = {0 , .01, .05, .10, .25, .50, .75, .90, .95, .99, 1};
Estimate = {80, 128, 151, 162, 184, 209, 236, 265, 284, 333, 727};
point=t(remove(Estimate,nrow(Quantile)));
width=t(remove(dif(Estimate),1));
prob=t(remove(dif(Quantile),1));
print Quantile Estimate point width prob;
/*Start to simulate uniform data*/
n=10000; /* sample size */
call randseed(123456789);
idx=t(sample(1:nrow(point),n,'replace',prob));
widths=width[idx];
points=point[idx];
eps=randfun(n,'uniform');
data=points+ceil(widths#eps); /* This is the simulate data */
/* Check the simulate data */
call histogram(data) label="Simulation" ;
call qntl(SimEst, data, Quantile);
print Quantile[F=percent6.] Estimate SimEst[F=5.];
quit;

 

View solution in original post

2 REPLIES 2
Ksharp
Super User

Check @Rick  's blog 

http://blogs.sas.com/content/iml/2014/06/18/distribution-from-quantiles.html 

 

 

Here is mine 

 

data have;
infile cards truncover expandtabs;
input Quantile Estimate ;
prob=dif(Quantile);
width=dif(Estimate);
cards;
0 80
.01 128
.05 151
.10 162
.25 184
.50 209
.75 236
.90 265
.95 284
.99 333
1 727
;
run;
proc sql;
select prob into : probs separated by ','
from have
where prob is not missing;
select Estimate into : points separated by ','
from have
where Quantile ne 1;
select width into : widths separated by ','
from have
where width is not missing;
select count(*)-1 into : count from have;
quit;
%let n=10000; /* sample size */
data want;
call streaminit(1234);
array points(&count) _temporary_ (&points);
array widths(&count) _temporary_ (&widths);
do i=1 to &n;
n=rand('table',&probs);
x=points{n}+ceil(widths{n}*rand('uniform'));
output;
end;
drop i;
run;
proc univariate data=want;
var x;
histogram x/kernel;
run;

 

 

 

 

Or An fast way is IML code , also written by me . Hope you like it .

 

proc iml;
/*
quantiles of total cholesterol from NHANES study
http://blogs.sas.com/content/iml/2014/06/18/distribution-from-quantiles.html
*/
Quantile = {0 , .01, .05, .10, .25, .50, .75, .90, .95, .99, 1};
Estimate = {80, 128, 151, 162, 184, 209, 236, 265, 284, 333, 727};
point=t(remove(Estimate,nrow(Quantile)));
width=t(remove(dif(Estimate),1));
prob=t(remove(dif(Quantile),1));
print Quantile Estimate point width prob;
/*Start to simulate uniform data*/
n=10000; /* sample size */
call randseed(123456789);
idx=t(sample(1:nrow(point),n,'replace',prob));
widths=width[idx];
points=point[idx];
eps=randfun(n,'uniform');
data=points+ceil(widths#eps); /* This is the simulate data */
/* Check the simulate data */
call histogram(data) label="Simulation" ;
call qntl(SimEst, data, Quantile);
print Quantile[F=percent6.] Estimate SimEst[F=5.];
quit;

 

DaLack
Calcite | Level 5

Many thanks for your support.

SAS Innovate 2025: Register Now

Registration is now open for SAS Innovate 2025 , our biggest and most exciting global event of the year! Join us in Orlando, FL, May 6-9.
Sign up by Dec. 31 to get the 2024 rate of just $495.
Register now!

What is Bayesian Analysis?

Learn the difference between classical and Bayesian statistical approaches and see a few PROC examples to perform Bayesian analysis in this video.

Find more tutorials on the SAS Users YouTube channel.

SAS Training: Just a Click Away

 Ready to level-up your skills? Choose your own adventure.

Browse our catalog!

Discussion stats
  • 2 replies
  • 1668 views
  • 0 likes
  • 2 in conversation