Fitting a Scaled Tweedie Model with Regressors
/*--------------------------------------------------------------
SAS Sample Library
Name: hsevex04.sas
Description: Example program from SAS/ETS High Performance
Procedures User's Guide, The HPSEVERITY Procedure
Title: Fitting a Scaled Tweedie Model with Regressors
Product: SAS High Performance Econometrics Software
Keys: Severity Distribution Modeling
PROC: HPSEVERITY
Notes:
--------------------------------------------------------------*/
/*--- Simulate a Tweedie sample that is affected by regressors ---*/
%let samplesize=300;
data test_sevtw(keep=y x1-x3
label='A Tweedie Sample Affected by Regressors');
array x{*} x1-x3;
array b{4} _TEMPORARY_ (5 0.25 -1 3);
call streaminit(45678);
label y='Response Influenced by Regressors';
Phi = 0.5;
P = 1.75;
do n = 1 to &samplesize;
Mu = 0;
do i = 1 to dim(x);
x(i) = rand('UNIFORM');
Mu = Mu + b(i+1) * x(i);
end;
Mu = b(1) * exp(Mu); /* b(1) is base value of mean */
cdf = rand('UNIFORM');
y = quantile('TWEEDIE', cdf, P, Mu, Phi);
output;
end;
run;
/*--- Fit the scale parameter version of the Tweedie distribution ---*/
proc hpseverity data=test_sevtw outest=estw covout print=all;
loss y;
scalemodel x1-x3;
dist stweedie;
nloptions tech=quanew;
run;
/*--- Compute the estimate, standard error, and p-value of the mean ---*/
data Mu0(keep=Parameter Estimate Stderr T probT);
Parameter='Mu0';
label Estimate='Estimate' Stderr='Standard;Error'
T='t Value' probT='Approx;Pr > |t|';
array cov{3,3} _temporary_;
array parms{3} _temporary_;
array parmrow{3} theta lambda p;
set estw(where=(_MODEL_='stweedie' and
(_type_='COV' or _type_='EST'))) end=eof;
n = _n_-1;
if (n = 0) then do;
do i=1 to 3;
parms[i] = parmrow[i];
end;
end;
else do;
do i=1 to n;
cov[n,i] = parmrow[i];
end;
end;
if (eof or n = 3) then do;
mu0 = parms[1]*parms[2]*(2-parms[3])/(parms[3]-1);
Estimate = mu0;
a = mu0/parms[1]; /* dMu0_T0 */
b = mu0/parms[2]; /* dMu0_L */
c = -mu0/((parms[3]-1)*(2-parms[3])); /* dMu0_P */
varMu0 = a * (a*cov[1,1]+b*cov[2,1]+c*cov[3,1]) +
b * (a*cov[2,1]+b*cov[2,2]+c*cov[3,2]) +
c * (a*cov[3,1]+b*cov[3,2]+c*cov[3,3]);
Stderr = sqrt(varMu0);
df = &samplesize-6;
T = mu0/Stderr;
probT = (1-probt(T, df))*2;
output;
stop;
end;
run;
proc print data=Mu0 noobs s=';';
run;