I have been searching how to create a bar graph with CI based on age adjusted proportions calculated through the command prop. The only related info that i found and was replicated in multiple forums was this explanation on means using the comman collapse https://stats.idre.ucla.edu/stata/fa...th-error-bars/. I think it does not apply for my case as I am calculating proportions. I solved the problem myself, but I think my solution is not very efficient, specially if you need to compare multiple groups. I'm putting the code I'm using below. Is there a simpler way to do this?
Code:
* Example generated by -dataex-. To install: ssc install dataex
clear
input float disease byte exp int agecat double ffinal_ajustaf
0 0 1 110.6365524
0 0 1 160.4600074
1 0 2 924.4567304
0 1 2 65.135256
0 0 1 108.1819536
0 0 1 121.648987
1 0 2 58.18526833
0 0 1 159.9585768
0 0 1 97.47475771
0 0 1 41.61796647
0 0 1 41.61796647
0 0 1 74.12985143
0 0 1 41.61796647
0 0 1 65.92337796
0 0 2 39.69179662
0 0 1 65.92337796
0 0 2 39.69179662
0 0 2 50.20095583
0 0 2 52.83390127
0 0 1 41.61796647
1 1 3 34.37288155
0 0 2 754.4944992
0 1 1 885.6033675
0 1 1 442.8016838
0 1 1 701.4033898
0 0 2 137.9927384
0 0 1 217.3974294
0 0 2 207.3358043
0 0 1 193.6139666
0 0 2 138.9096722
0 0 2 155.5018533
0 0 1 258.2701797
0 0 1 190.9400274
0 0 1 2469.478651
0 0 2 250.870378
0 0 2 313.6776802
0 0 2 396.729821
0 0 1 200.7059928
0 1 2 305.6770752
0 0 1 126.7073311
0 0 1 376.3553746
0 0 2 686.8426124
0 0 1 376.3553746
0 0 2 447.6729288
0 0 1 524.6839447
0 0 2 102.8660223
0 0 1 86.4787166
0 0 1 271.1393564
0 1 2 140.9825176
0 0 2 216.5223021
0 1 1 172.9574332
0 0 1 241.1231363
0 0 1 138.9886564
0 0 2 118.412125
0 0 2 176.4457413
0 0 1 162.7648676
0 0 1 277.9773127
0 0 2 441.0522683
0 0 2 419.0726959
0 0 1 347.4227357
0 0 2 882.1045366
0 0 1 347.4227357
1 0 3 286.9414716
0 0 2 832.4822473
0 0 1 808.2117563
0 0 2 383.3131622
0 0 2 257.2401337
0 0 3 334.3509851
0 0 1 301.9408742
0 0 1 119.9689326
0 0 2 87.2779025
0 1 2 97.70288399
1 0 4 113.4405128
0 0 1 2108.744502
0 0 1 3340.277591
0 1 3 1741.64264
0 0 1 3340.277591
0 0 1 7512.184273
0 0 2 2011.147228
0 0 1 434.9386402
0 0 1 774.7119696
0 0 1 169.7805946
1 0 4 119.7407084
0 0 1 144.9795467
0 0 1 774.7119696
0 0 2 77.13385578
0 0 3 89.55829958
0 0 1 189.4246304
0 0 1 24.41371755
0 0 1 14.95384761
0 0 1 35.02388022
0 0 1 48.82743511
0 0 3 24.41369597
0 0 3 24.41369597
1 1 3 24.41369597
1 0 3 24.41369597
1 0 3 24.41369597
0 0 1 256.8417883
0 0 3 123.7009242
0 0 1 201.3209721
end
label values disease disease
label def disease 0 "Normal", modify
label def disease 1 "Disease", modify
label values exp exp
label def exp 0 "No", modify
label def exp 1 "Yes", modify
label values agecat agecat
label def agecat 1 "60-69", modify
label def agecat 2 "70-79", modify
label def agecat 3 "80-89", modify
label def agecat 4 "90-99", modify
label values disease disease
label def disease 0 "Normal", modify
label def disease 1 "Disease", modify
label values exp exp
label def exp 0 "No", modify
label def exp 1 "Yes", modify
label values agecat agecat
label def agecat 1 "60-69", modify
label def agecat 2 "70-79", modify
label def agecat 3 "80-89", modify
svyset [pw=ffinal_ajustaf]
mat def cdpn = J(7, 1, 0)
mat def cdpa = J(7, 1, 0)
mat def cdp_ub = J(7, 1, 0)
mat def cdp_lb = J(7, 1, 0)
svy: tab agecat
mat def b= e(b)
mat def b= b'
gen agestd_wt=.
replace agestd_wt= b[1,1] if agecat==1
replace agestd_wt= b[2,1] if agecat==2
replace agestd_wt= b[3,1] if agecat==3
replace agestd_wt= b[4,1] if agecat==4
replace agestd_wt= b[5,1] if agecat==5
local i=2
local j=1
forval n=0/1{
svy: prop disease if exp==`n', stdize(agecat) stdweight(agestd_wt) over (exp)
matrix def r= r(table)
matrix def r=r'
matrix cdpn[`i',1]= r[`j',1]
mat cdp_ub[`i',1] = r[`j',6]
mat cdp_lb[`i',1] = r[`j',5]
local i=`i'+1
local j=`j'+1
matrix cdpa[`i',1]= r[`j',1]
mat cdp_ub[`i',1] = r[`j',6]
mat cdp_lb[`i',1] = r[`j',5]
local i=`i'+2
local j=1
}
mat table= cdpn, cdpa, cdp_lb, cdp_ub
mat colname table= Prevalence lb ub
svmat table
gen i=_n
gen ibar=i in 1/7
replace table1=. if table1==0
replace table2=. if table2==0
replace table3=. if table3==0
replace table4=. if table4==0
twoway (bar table1 table2 ibar) || rcap table3 table4 ibar, xlabel(2.5 "No exposure" 5.5 "Exposure") xtitle("")ytitle("Proportion") ylabel(0(.2)1) legend(order(1 "Normal" 2 "Disease")) title("Age Adjusted Prevalence of Disease")
