/*Fasting Weight Adjusted*/

/*Program highlights:
1. First uses the Fasting weight because fasting plasma glucose measurements are only available
for the fasting population who provided blood at the clinic.
2. An Adjustment Weight for the Fasting Weight was then created due to missing data in the
subsetted sample.
3. Calculates Crude prevalence estimates of diabetes and corresponding Relative 
Standard Errors (RSEs).
4. Calculates Age-Adjusted prevalence estimates of diabetes and corresponding RSEs.
*/

proc format;

/*gender*/
value riagendr
	1='Male'
	2='Female';

/*Race_eth*/
value race_eth
	1='Non-Hisp White'
	2='Non-Hisp Black'
	3='Non-Hisp Asian' /*(Non-Hispanic Asian and Hawaiian/PI)*/
	4='Hispanic'
	5='Non-Hisp Other' /*(American Indian, Other, Multiracial who could not select one race, and
	SPs missing race but Non-Hispanic)*/;

/*Total Diabetes*/
	value Dbts
		1='Diabetes' /*FG >=126 or Previous Dx*/
		2='Non-Diabetes'; /*FG <126 and Never Told*/

/*Self-Report Prev Diagnosed*/
	value SRDbts
		1='Previous Dx'
		2='Never told';

/*Undiagnosed Among Total Population*/
	value DxDbts
		1='Previous Dx'
		2='Undiagnosed'
		3='Non-Diabetic';

/*Undiagnosed Among Diabetics*/
	value DxDbts2f
		1='Previous Dx'
		2='Undiagnosed';

run;


libname hanes 'C:\Data';

proc sort data=hanes.spfile out=spfile; by sp_id; run;
proc sort data=hanes.capi out=capi; by sp_id; run;
proc sort data=hanes.labs out=labs; by sp_id; run;

data diabetes; merge spfile (in=a) capi labs; by sp_id;

if a and wtsf1f>0 and diq010 ne .; /*Limits complete SPs in Fasting Sample with diabetes 
	interview data*/
	/*Do not need to specify limitations for missing lbxglu - in our dataset, only SPs who are
	self-reported diabetics are missing lbxglu (N=4) and are classified as diabetic regardless
	of lbxglu value.*/

/*Total Diabetes*/
	if lbxglu ge 126 or diq010=1 then Dbts=1; /*Fasting Glucose >= 126 or Previously Diagnosed*/
	else if lbxglu lt 126 and diq010 ne 1 then Dbts=2; /*Fasting Glucose <126 and Never told/
		Refused/Missing in DIQ010*/

/*Undiagnosed among Total Population*/
	if Dbts=1 and diq010=1 then DxDbts=1; /*Previously Diagnosed*/
	else if Dbts=1 and diq010 ne 1 then DxDbts=2; /*Undiagnosed*/
	else if Dbts ne 1 and diq010 ne 1 then dxdbts=3; /*Non-Diabetic*/

/*Undiagnosed among Diabetics*/
	if Dbts=1 and diq010=1 then DxDbts2=1; /*Previously Diagnosed*/
	else if Dbts=1 and diq010 ne 1 then DxDbts2=2; /*Undiagnosed*/

format dbts dbts. riagendr riagendr. race_eth race_eth. dxdbts dxdbts. dxdbts2 dxdbts2f.;

run;


/*Merge in Adjustment Weights for 'FASTING' weight - Refer to Analytic Guidelines for
	further details.*/

/*Sum Fasting weights of the respondents (in subsetted dataset) within each category of age,
gender, race/ethnicity - 'SUM B'*/
proc summary data=diabetes nway;
var wtsf1f; /*Fasting Weight*/
class agewt riagendr racewt;
output out=samplewts sum=samplewt;
run;

/*The file hanes.fst_con_tots contains the sum of the Fasting Control Totals within each
category of age, gender, race/ethnicity - 'SUM A'*/
/*Merge SUM B to SUM A to create Weight Factor (A/B)*/
proc sort data=samplewts; by agewt riagendr racewt;
proc sort data=hanes.fst_con_tots out=fst_con_tots; by agewt riagendr racewt;
data wtfactor; merge samplewts fst_con_tots; by agewt riagendr racewt;
wtfactor=fasting_control_total/samplewt; /* Weight Factor = A/B */
run;

/*Merge dataset with new weight factor to subsetted dataset in order to adjust the weight by the
weight factor A/B for each SP in the cell*/
proc sort data=wtfactor; by agewt riagendr racewt;
proc sort data=diabetes; by agewt riagendr racewt; run;
data diabad; merge diabetes wtfactor; by agewt riagendr racewt;
diabwt = wtsf1f * wtfactor; /*New weight = Old weight * (A/B)*/
run;
proc sort data=diabad; by stratum psu; run;

/*Check adjusted weight*/
proc summary data=diabad nway;
var diabwt;
output out=popchk sum=pop;
run;
proc print; run; /*5,827,719*/


/*Prevalence of Diabetes - CRUDE (using PROC CROSSTAB)*/
proc sort data=diabad; by stratum psu; run;
PROC CROSSTAB data=diabad;
NEST STRATUM PSU;
WEIGHT diabwt; /* Adjusted Fasting Weight */
subgroup dbts riagendr agegroup race_eth;
LEVELS     2     2     	 5		  5;
TABLES  RIAGENDR*dbts;
RTITLE "Diabetes Prevalence Among NYC Population, by Gender - Crude";
PRINT NSUM="SAMSIZE" WSUM="POPSIZE" ROWPER SEROW lowrow uprow/
TESTS=ALL WSUMFMT=F9.0 SEROWFMT= F6.3 style=nchs;
output/filename=dbcrude filetype=sas TABLECELL=DEFAULT REPLACE;
SETENV COLWIDTH=14 DECWIDTH=2 COLSPCE=2; run;

/*Calculate Relative Standard Error of the estimates using the SUDAAN output dataset.  Typically,
estimates with RSE >=30% are considered unreliable and notated that it should be interpreted with
caution.*/
data rse; set dbcrude;
RSE=round(serow/rowper*100,0.1);
run;


/*Prevalence of Diabetes - Age-Adjusted (using PROC DESCRIPT)*/
proc sort data=diabad; by stratum psu; run;
PROC DESCRIPT data=diabad design=wr;
nest stratum psu;
weight diabwt; /* Adjusted Fasting Weight */
VAR dbts;
catlevel 1;
subgroup ageadj _one_;
LEVELS      3	  1;
tables _one_;
stdvar ageadj;
stdwgt 0.396579 0.371795 0.231626;
print nsum wsum="weights for totals" total="weights for %" percent sepercent lowpct uppct/
style=nchs wsumfmt=F8.0 totalfmt=F8.0 percentfmt=f6.3 lowpctfmt=f6.3 uppctfmt=f6.3 sepercentfmt=f6.3;
output/filename=dbagead filetype=sas TABLECELL=DEFAULT percentfmt=f8.3 lowpctfmt=f8.3 uppctfmt=f8.3
sepercentfmt=f6.3 REPLACE;
RTITLE "Diabetes NYC - Age-Adjusted";
SETENV COLWIDTH=16 DECWIDTH=2 COLSPCE=2;
run;

/*Calculate Relative Standard Error of the estimates using the SUDAAN output dataset.  Typically,
estimates with RSE >=30% are considered unreliable and notated that it should be interpreted with
caution.*/
data DiabTabs; set dbagead;
RSE=round(sepercent/percent*100,.001);
keep _one_ variable NSUM WSUM total percent sepercent lowpct uppct rse;
run;