****pgm: linear_stata.do******** log using "c:\NHANES\log\linear.log", replace use "C:\NHANES\DATA\analysis_data.dta", clear ***create variable codes***** capture label drop educ gen educ= dmdeduc replace educ=. if educ >3 gen smoker=1 if smq020==2 replace smoker=2 if smq020==1 & smq040==3 replace smoker=3 if smq020==1 & (smq040==1 | smq040==2) gen bmicat=1 if bmxbmi>=0 & bmxbmi<18.5 replace bmicat=2 if bmxbmi>=18.5 & bmxbmi<25 replace bmicat=3 if bmxbmi>=25 & bmxbmi<30 replace bmicat=4 if bmxbmi>=30 & bmxbmi<. /*eligible is 1 if other variables used in final regression model are non-blank; this keeps #obs same*/ gen eligible=1 if (ridageyr >=20 & ridageyr <.) & (bmxbmi !=. & riagendr !=. & ridreth1 !=. & ridageyr !=. & smoker !=. & educ !=.) tab bmicat if eligible==1, gen(Ibmicat) ****format variables****** capture label drop sexfmt label define sexfmt 1 "male" label define sexfmt 2 "female", add label define race2fmt 1 "Mex American" label define race2fmt 3 "NH White", add label define race2fmt 4 "NH Black", add capture label drop smkfmt label define smkfmt 1 "Never smoker" label define smkfmt 2 "Past smoker", add label define smkfmt 3 "Current smoker", add label define educ 1 "< High school" label define educ 2 "High school", add label define educ 3 "> High school", add label define bmicat 1 "underweight" label define bmicat 2 "normal weight", add label define bmicat 3 "overweight", add label define bmicat 4 "obese", add label values riagendr sexfmt label values ridreth1 race2fmt label values smoker smkfmt label values educ educ label values bmicat bmicat ****specify survey design variables**** svyset sdmvpsu [pweight=wtmec4yr], strata(sdmvstra) vce(linearized) ***simple linear regression*** svy, subpop(if eligible==1): regress lbdhdl bmxbmi svy, subpop(if eligible==1): regress lbdhdl bmicat svy, subpop(if eligible==1): mean lbdhdl, over(bmicat) svy, subpop(if eligible==1): mean lbdhdl, over(riagendr) svy, subpop(if eligible==1): mean lbdhdl, over(ridreth1) svy, subpop(if eligible==1): mean lbdhdl, over(smoker) svy, subpop(if eligible==1): mean lbdhdl, over(educ) ****change reference level for categorical variables in model****** char ridreth1[omit]3 char smoker[omit]3 char educ[omit]3 char bmicat[omit]2 ****NOTE: xi and i.var are used to expand and denote categorical variables***** xi:svy, subpop(if eligible==1): regress lbdhdl i.bmicat ****multiple linear regression*** xi: svy, subpop(if ridageyr >=20 & ridageyr <.): regress lbdhdl /// i.riagendr i.ridreth1 ridageyr i.smoker i.educ i.bmicat *** NOTE: The -adjust- command uses only the sample mean, not the mean based on the survey design, *** *** when performing its computations. Therefore, if you want to use the survey mean, *** *** you would need to calculate it first and specify it explicitly in the -adjust- command. *** *** The following commands use -summarize- which is an rclass command and will not cause any *** *** trouble if run between the -svy: regress- and -adjust- commands whereas -svy:mean- is an *** *** eclass command and cannot be used in between these commands. *** ****get survey means within the subpopulation**** quietly { sum _Iridreth1_1 [aw=wtmec4yr] if ridageyr >=20 & ridageyr <. & e(sample) local rid1 = r(mean) sum _Iridreth1_2 [aw=wtmec4yr] if ridageyr >=20 & ridageyr <. & e(sample) local rid2 = r(mean) sum _Iridreth1_4 [aw=wtmec4yr] if ridageyr >=20 & ridageyr <. & e(sample) local rid4 = r(mean) sum _Iridreth1_5 [aw=wtmec4yr] if ridageyr >=20 & ridageyr <. & e(sample) local rid5 = r(mean) sum ridageyr [aw=wtmec4yr] if ridageyr >=20 & ridageyr <. & e(sample) local ridage = r(mean) sum _Ismoker_1 [aw=wtmec4yr] if ridageyr >=20 & ridageyr <. & e(sample) local smoke1 = r(mean) sum _Ismoker_2 [aw=wtmec4yr] if ridageyr >=20 & ridageyr <. & e(sample) local smoke2 = r(mean) sum _Ieduc_1 [aw=wtmec4yr] if ridageyr >=20 & ridageyr <. & e(sample) local educ1 = r(mean) sum _Ieduc_2 [aw=wtmec4yr] if ridageyr >=20 & ridageyr <. & e(sample) local educ2 = r(mean) sum _Ibmicat_1 [aw=wtmec4yr] if ridageyr >=20 & ridageyr <. & e(sample) local bmicat1 = r(mean) sum _Ibmicat_3 [aw=wtmec4yr] if ridageyr >=20 & ridageyr <. & e(sample) local bmicat3 = r(mean) sum _Ibmicat_4 [aw=wtmec4yr] if ridageyr >=20 & ridageyr <. & e(sample) local bmicat4 = r(mean) sum _Iriagendr_2 [aw=wtmec4yr] if ridageyr >=20 & ridageyr <. & e(sample) local riagendr2 = r(mean) } ****Use adjust command to produced adjusted means**** adjust _Iridreth1_1=`rid1' _Iridreth1_2=`rid2' _Iridreth1_4=`rid4' /// _Iridreth1_5=`rid5' ridageyr=`ridage' _Ismoker_1=`smoke1' /// _Ismoker_2=`smoke2' _Ieduc_1=`educ1' _Ieduc_2=`educ2' /// _Ibmicat_1=`bmicat1' _Ibmicat_3=`bmicat3' _Ibmicat_4=`bmicat4' /// if ridageyr >=20 & ridageyr <. & e(sample), by(riagendr) se adjust _Iriagendr_2=`riagendr2' ridageyr=`ridage' _Ismoker_1=`smoke1' /// _Ismoker_2=`smoke2' _Ieduc_1=`educ1' _Ieduc_2=`educ2' /// _Ibmicat_1=`bmicat1' _Ibmicat_3=`bmicat3' _Ibmicat_4=`bmicat4' /// if ridageyr >=20 & ridageyr <. & e(sample), by(ridreth1) se adjust _Iriagendr_2=`riagendr2' _Iridreth1_1=`rid1' _Iridreth1_2=`rid2' /// _Iridreth1_4=`rid4' _Iridreth1_5=`rid5' ridageyr=`ridage' /// _Ieduc_1=`educ1' _Ieduc_2=`educ2' _Ibmicat_1=`bmicat1' /// _Ibmicat_3=`bmicat3' _Ibmicat_4=`bmicat4' if ridageyr >=20 & /// ridageyr <. & e(sample), by(smoker) se adjust _Iriagendr_2=`riagendr2' _Iridreth1_1=`rid1' _Iridreth1_2=`rid2' /// _Iridreth1_4=`rid4' _Iridreth1_5=`rid5' ridageyr=`ridage' /// _Ismoker_1=`smoke1' _Ismoker_2=`smoke2' _Ibmicat_1=`bmicat1' /// _Ibmicat_3=`bmicat3' _Ibmicat_4=`bmicat4' if ridageyr >=20 & /// ridageyr <. & e(sample), by(educ) se adjust _Iriagendr_2=`riagendr2' _Iridreth1_1=`rid1' _Iridreth1_2=`rid2' /// _Iridreth1_4=`rid4' _Iridreth1_5=`rid5' ridageyr=`ridage' /// _Ieduc_1=`educ1' _Ieduc_2=`educ2' _Ismoker_1=`smoke1' /// _Ismoker_2=`smoke2' if ridageyr >=20 & ridageyr <. & e(sample), by(bmicat) se ***Model using an interaction term*** **xi: svy, subpop(if ridageyr >=20 & ridageyr <.): regress lbdhdl i.riagendr i.ridreth1 ridageyr i.smoker i.educ i.smoker*bmxbmi** ****test all coeffiecients in model together**** **test** ****test each coefficient in model separately*** **test _Iriagendr_2, nosvyadjust** **test _Iridreth1_1 _Iridreth1_2 _Iridreth1_4 _Iridreth1_5, nosvyadjust** **test ridageyr, nosvyadjust** **test _Ibmicat_1 _Ibmicat_3 _Ibmicat_4, nosvyadjust** **test _Ismoker_1 _Ismoker_2, nosvyadjust** **test _Ieduc_1 _Ieduc_2, nosvyadjust** *****test never smokers vs. past smokers**** **test _Ismoker_1 - _Ismoker_2 =0** ***Or, can use: test _Ismoker_1 = _Ismoker_2 ****** ********* **log close**