Html Version of SAS Code to run exemplar 1
Analyses from total scores.
It is intended to show you the code and to allow links, not to use
as a SAS program. The SAS program is the file ex1.sas which you should
save to a file and/or read in to SAS.
To see output from the commands go to the SAS results file.
Links in this page
Mean income with different design assumptions
Subgroup lone parents
Percentiles
Raking with CALMAR
Jacknife macro for mean
EVERYTHING INSIDE /* AND */ IS TAKEN AS A COMMENT in SAS programs
These are shown in green here
back to top 
/****************************************************************************************;
ex1_hbai
analysis for peas project exemplar 1
Gillian Raab 5/04
data ex1.ex1; set ex1.ex1; drop ad12 fyear;run;
*********************************************************************/
libname ex1 'C:\Documents and Settings\gillian raab\My Documents\aprojects\peas\ex1datafiles\exemp1\data';
/*---------------------------------------------
first analyses of mean income per household
-----------------------------------------------*/
title 'Mean income without weighting from PROC MEANS';
proc means data=ex1.ex1 mean stderr; * inappropriate assumption of total random sampling;
var hhinc; run;
title 'Mean income with weighting from PROC MEANS';
proc means data=ex1.ex1 mean stderr; * inappropriate assumption of total random sampling;
var hhinc;
weight gross2; run;
title 'Survey means analysis with weighting from PROC SURVEYMEANS';
proc surveymeans data=ex1.ex1 mean stderr;
weight gross2;
var hhinc;
run;
title 'Mean income with weighting and clustering from PROC SURVEYMEANS';
proc surveymeans data=ex1.ex1 mean stderr; * introducing clusters ;
weight gross2;
cluster psu;
var hhinc;
run;
/*------------------------------------------------------------
note that this next analysis is wrong as it makes too many psus
by splitting up the clusters by strata
-------------------------------------------------------------*/
title 'Survey means with post-strata wrongly defined as strata';
proc surveymeans data=ex1.ex1 mean stderr; * introducing clusters and strata by ctband ;
weight gross2;
cluster psu;
var hhinc ;
strata ctband;
run;
/*--------------------------------------------------------------
now return to previous model without post strata
and get mean income by subgroup
-----------------------------------------------------------------*/
data ex1new; set ex1.ex1;
if adulth=1 and depchldh>0 then lonep=1;
else lonep=0;
run;
title 'income for lone parents';
proc freq data=ex1new; table lonep; run;
proc surveymeans data=ex1new mean stderr; * introducing clusters ;
weight gross2;
cluster psu;
domain lonep;
var hhinc;
run;
/*-----------------------------------------------------------------
now save the file containing the macro (pctilegrps.sas) from the
web page to your computer and change the code below so that it
includes the file from YOUR machine
You can open it to get help in the comments field on how to use
------------------------------------------------------------------*/
back to top
%include "C:\Documents and Settings\gillian raab\My Documents\aprojects\peas\ex1datafiles\exemp1\program_code\pctilegrps.sas";
/*--------------------------------------------------------------------
now use it to get the s.e.s and c.i.s for percentiles in one group.
These will be somewhat too wide as they do not allow for the
precision gained by post-stratification. This seems to make more difference
at the middle percentiles.
-------------------------------------------------------------*/
proc format ; value lonep 0='not lone parent' 1='lone parent';run;
/*---------------- first percentiles all-----------------------*/
%pctilegrps(var=hhinc,
data=ex1new,
group=one,
n=1,
weight=gross2,
strata=,
cluster=psu,
final= new);
/*----------------and now grouped by lone parents----------------------*/
%pctilegrps(var=hhinc,
data=ex1new,
group=lonep,
n=2,
weight=gross2,
strata=,
format =format lonep lonep.,
cluster=psu,
final= new);
back to top
/*--------------------------------------------------------------------
Next code rakes survey using the MACRO CALMAR
available from INSEE
http://www.insee.fr/fr/nom_def_met/outils_stat/calmar/cal_res.htm
Documentation is in French. Key words to know
POIDS = weights
OUI= yes
NON = no
To use download tyhe ZIP file and store the file called
SASMACR.SAS7Bcat it in a SAS data library.
Here it is in the same directory as the data from exemplar 6.
Set options to allow SAS to se it
---------------------------------------------------------------------*/
options nomprint;
options mstored sasmstore=ex1;
/*------------------------------------------------------
now set up table of margins
------------------------------------------------------------*/
data margins;
input var $ n mar1 mar2 mar3 mar4 mar5 mar6 mar7 mar8 mar9;
CARDS;
CTBAND 9 24.83 24.62 15.45 11.91 11.96 5.95 3.94 0.45 0.89
tenure 4 62.63 21.59 5.58 10.20 . . . . .
;
/*----------------------------------------------------------
now rake data with CALMAR M=2 gives the method of raking
obseli=OUI (yes) sends observations eliminated to another file
effpop is the effective population size
pct=OUI indicates data in %s
--------------------------------------------------------------*/
TITLE 'RAKING WITH CALMAR';
%calmar(data=ex1.ex1,poids=GROSS2,ident=SERNUM,datamar=margins,M=2,
editpoi=oui,obseli=oui,datapoi=newdata,poidsfin=newweight,
labelpoi='new weight',pct=oui,effpop=2242012);
proc sort data=ex1.ex1; by sernum;
proc sort data=newdata; by sernum;* just new weight and ID;
data poststrat; * this is the new complete data set with the raked weights;
merge ex1.ex1 newdata;
by sernum;
run;
proc format ;
value tenure 1='owned' 2='rented LA' 3='oth pub rent' 4='private rent or rent free';
value CTBAND 1='A' 2='B' 3='C' 4='D' 5='E' 6='F' 7='G'8='H' 9='not sep valued';
run;
proc gplot data=poststrat; * check to see how much weights have changed;
plot newweight*gross2=tenure;
plot newweight*gross2=ctband;
format tenure tenure. ;
run;
proc surveymeans data=poststrat mean stderr; * now get new survey mean analysis ;
weight newweight;
cluster psu;
var hhinc;
run;
back to top
/*--------------------------------------------------------
this next bit of code does a Jacknife analysis
It is not meant to be comprehensive, just to give an
indication of how it could be done. It is not
efficient.
Send output to a file and switch off listing before you run it
so as not to clog system
----------------------------------------------------*/
ods html file='temp.htm'
path='c:\' (url=none)
style=minimal;
/*------------------- first make the jacknife weights and rake each one of them-------------*/
%macro makejkwts(data=ex1.ex1,origwt=GROSS2,out=jck,psu=PSU,npsu=320,
effpop=2242012,sernum=SERNUM,M=2,margins=margins);
* first rake the original data set;
%calmar(data=&data,poids=&origwt,ident=&sernum,datamar=&margins,M=&M,
editpoi=non,obseli=oui,datapoi=alldata,poidsfin=newweight,
pct=oui,effpop=&effpop,stat=non);
data &out; merge alldata &data; run;
proc sort data=&out; by sernum;
%do nrep=1 %to &npsu;
data new; set &data;
if &psu=&nrep then &origwt=0;
run;
%calmar(data=new,poids=&origwt,ident=&sernum,datamar=&margins,M=&M,
editpoi=non,obseli=oui,datapoi=newdata,poidsfin=newweight&nrep,
pct=oui,effpop=&effpop,stat=non);
data new;
set newdata __obseli (in=excl);
if excl then newweight&nrep=0;
keep &sernum newweight&nrep;
run;
proc sort data=new; by &sernum; run;
data &out ;merge &out new; by sernum; run;
%end;
%mend;
%makejkwts(out=ex1.ex1repwts);
ods html close;
/*----------------- now a macro to do a simple jacknife------------------*/
/******************************************
input parameters:
indata - input data set.
Outdata - output data (mean and std err.)
var - analysis variable on input data set.
wt - variable name on input data for full sample weights.
Repwt - prefix used in replicate weight names. Replicate weight
names must be in the form: repwt01, repwt02, ... repwt45
n - number of replicate weights on input data set.
Program creates an output dataset containing the method, the
number of replicates, the mean, and the standard error.
********************************************************/
%macro jack2_mean(indata,outdata,var,repwt,n);
%do i=1 %to &n; ** cycle through the replicate weights **;
%put calculating results for replicate weight &i;
proc means data=&indata noprint;
var &var;
weight &repwt&i;
output out=rep_mean_junk mean=mean;
run;
data rep_mean_junk2; set rep_mean_junk2
rep_mean_junk(in=j);
if j then rep=&i;
run;
%end;
proc means data=rep_mean_junk2 noprint;
var mean;
output out=&outdata mean=mean var=var;
run;
data &outdata; set &outdata(keep=mean var);
_METHOD_ = "JK2";
_REPLICATES_=&n;
&var._mean = mean;
&var._std = sqrt(var*(&n-1)*(&n-1)/&n);
*drop mean var;
run;
proc print;
run;
%mend jack2_mean;
%jack2_mean(ex1.ex1repwts,res,hhinc,newweight,320);