## Defining a Mixed-Tail Distribution

```/*--------------------------------------------------------------

SAS Sample Library

Name: hsevex03.sas
Description: Example program from SAS/ETS High Performance
Procedures User's Guide, The HPSEVERITY Procedure
Title: Defining a Mixed-Tail Distribution
Product: SAS High Performance Econometrics Software
Keys: Severity Distribution Modeling
PROC: HPSEVERITY
Notes:

--------------------------------------------------------------*/

/*------- Define Lognormal Body-GPD Tail Mixed Distribution -------*/
proc fcmp library=sashelp.svrtdist outlib=work.sevexmpl.models;
function LOGNGPD_DESCRIPTION() \$256;
length desc \$256;
desc1 = "Lognormal Body-GPD Tail Distribution.";
desc2 = " Mu, Sigma, and Xi are free parameters.";
desc3 = " Xr and Pn are constant parameters.";
desc = desc1 || desc2 || desc3;
return(desc);
endsub;

function LOGNGPD_SCALETRANSFORM() \$3;
length xform \$3;
xform = "LOG";
return (xform);
endsub;

subroutine LOGNGPD_CONSTANTPARM(Xr,Pn);
endsub;

function LOGNGPD_PDF(x, Mu,Sigma,Xi,Xr,Pn);
cutoff = exp(Mu) * Xr;
p = CDF('LOGN',cutoff, Mu, Sigma);
if (x < cutoff + constant('MACEPS')) then do;
return ((Pn/p)*PDF('LOGN', x, Mu, Sigma));
end;
else do;
gpd_scale = p*((1-Pn)/Pn)/PDF('LOGN', cutoff, Mu, Sigma);
h = (1+Xi*(x-cutoff)/gpd_scale)**(-1-(1/Xi))/gpd_scale;
return ((1-Pn)*h);
end;
endsub;

function LOGNGPD_CDF(x, Mu,Sigma,Xi,Xr,Pn);
cutoff = exp(Mu) * Xr;
p = CDF('LOGN',cutoff, Mu, Sigma);
if (x < cutoff + constant('MACEPS')) then do;
return ((Pn/p)*CDF('LOGN', x, Mu, Sigma));
end;
else do;
gpd_scale = p*((1-Pn)/Pn)/PDF('LOGN', cutoff, Mu, Sigma);
H = 1 - (1 + Xi*((x-cutoff)/gpd_scale))**(-1/Xi);
return (Pn + (1-Pn)*H);
end;
endsub;

subroutine LOGNGPD_PARMINIT(dim,x[*],nx[*],F[*],Ftype,
Mu,Sigma,Xi,Xr,Pn);
outargs Mu,Sigma,Xi,Xr,Pn;
array xe[1] / nosymbols;
array nxe[1] / nosymbols;

eps = constant('MACEPS');

Pn = 0.8; /* Set mixing probability */
_status_ = .;
call streaminit(56789);
Xb = svrtutil_hillcutoff(dim, x, 100, 25, _status_);
if (missing(_status_) or _status_ = 1) then
Xb = svrtutil_percentile(Pn, dim, x, F, Ftype);

/* Initialize lognormal parameters */
call logn_parminit(dim, x, nx, F, Ftype, Mu, Sigma);
if (not(missing(Mu))) then
Xr = Xb/exp(Mu);
else
Xr = .;

/* prepare arrays for excess values */
i = 1;
do while (i <= dim and x[i] < Xb+eps);
i = i + 1;
end;
dime = dim-i+1;
if (dime > 0) then do;
call dynamic_array(xe, dime);
call dynamic_array(nxe, dime);
j = 1;
do while(i <= dim);
xe[j] = x[i] - Xb;
nxe[j] = nx[i];
i = i + 1;
j = j + 1;
end;

/* Initialize GPD's shape parameter using excess values */
call gpd_parminit(dime, xe, nxe, F, Ftype, theta_gpd, Xi);
end;
else do;
Xi = .;
end;
endsub;

subroutine LOGNGPD_LOWERBOUNDS(Mu,Sigma,Xi,Xr,Pn);
outargs Mu,Sigma,Xi,Xr,Pn;

Mu   = .; /* Mu has no lower bound */
Sigma = 0; /* Sigma > 0 */
Xi   = 0; /* Xi > 0 */
endsub;
quit;

/*----- Simulate a sample for the mixed-tail distribution -----*/
data testmixdist(keep=y label='Lognormal Body-GPD Tail Sample');
call streaminit(45678);
label y='Response Variable';
N = 100;
Mu = 1.5;
Sigma = 0.25;
Xi = 1.5;
Pn = 0.8;

/* Generate data comprising the lognormal body */
Nbody = N*Pn;
do i=1 to Nbody;
y = exp(Mu) * rand('LOGNORMAL')**Sigma;
output;
end;

/* Generate data comprising the GPD tail */
cutoff = quantile('LOGNORMAL', Pn, Mu, Sigma);
gpd_scale = (1-Pn) / pdf('LOGNORMAL', cutoff, Mu, Sigma);
do i=Nbody+1 to N;
y = cutoff + ((1-rand('UNIFORM'))**(-Xi) - 1)*gpd_scale/Xi;
output;
end;
run;

/*--- Set the search path for functions defined with PROC FCMP ---*/
options cmplib=(work.sevexmpl);
/*-------- Fit LOGNGPD model with PROC HPSEVERITY --------*/
proc hpseverity data=testmixdist print=all outest=parmest;
loss y;
dist logngpd burr logn gpd;
run;

/*-------- Compute tail cutoff and tail distribution's scale --------*/
data xb_thetat(keep=x_b theta_t);
set parmest(where=(_MODEL_='logngpd' and _TYPE_='EST'));
x_b = exp(Mu) * Xr;
theta_t = (CDF('LOGN',x_b,Mu,Sigma)/PDF('LOGN',x_b,Mu,Sigma)) *
((1-Pn)/Pn);
run;

proc print data=xb_thetat noobs;
run;

```