I am having problems extracting coefficients from xtgee model. What I am interested in extracting are _cons 1.p203 1.period 1.p203#1.period. 1.p203#1.period will be from the margins, contrast. Your assistance will be much appreciated. data is here attached.sample.dta
Code:
xtset participantid period
* Analysis section - corrected
local did_bin pp110_1 // pp110_2 pp211_1 pp211_2
local cov age edu mdp par hfd
foreach outcome of local did_bin {
* Clear existing collections
collect clear
di _n "Processing outcome: `outcome'"
* Run GEE model with verbose output
capture noisily xtgee `outcome' i.p203##i.period i.age i.edu i.mdp i.par i.hfd [pweight=psweight], ///
family(binomial) link(logit) corr(exchangeable) vce(robust) nolog level(95)
if _rc {
di as error "GEE failed for `outcome' (rc = " _rc "), skipping"
continue
}
* Check model convergence
di "Converged: " e(converged)
* Store coefficients and variance
matrix beta = e(b)
matrix V = e(V)
* Debug: Display matrix contents
di _n "Coefficient matrix (beta):"
matrix list beta
di _n "Variance matrix (V):"
matrix list V
* Create collection
collect create `outcome'_results
di "Collection `outcome'_results created"
* Define coefficient labels and variables
local coef_labels "Intercept p203 period interaction"
local coef_vars "_cons 1.p203 1.period 1.p203#1.period"
local i = 1
* Debug: Verify coef_vars and i
di "coef_vars: `coef_vars'"
di "Starting i: `i'"
* Add coefficients and CIs to collection
foreach coef in `coef_labels' {
local var : word `i' of `coef_vars'
di _n "Processing coefficient: `coef' (variable: `var', i = `i')"
* Check if column exists in beta and V
local col_idx = colnumb(beta, "`var'")
if `col_idx' == . {
di as error "Column `var' not found in beta matrix, skipping"
continue
}
local col_idx_V = colnumb(V, "`var'")
if `col_idx_V' == . {
di as error "Column `var' not found in V matrix, skipping"
continue
}
di "Column index (beta): `col_idx'"
di "Column index (V): `col_idx_V'"
* Extract estimate
capture scalar est = beta[1, `col_idx']
if _rc {
di as error "Failed to extract est for `var' (rc = " _rc "), skipping"
continue
}
if missing(est) {
di as error "est is missing for `var', skipping"
continue
}
* Extract variance and check for validity
capture scalar var_val = V[`col_idx_V', `col_idx_V']
if _rc {
di as error "Failed to extract variance for `var' (rc = " _rc "), skipping"
continue
}
if missing(var_val) | var_val < 0 {
di as error "Invalid variance for `var' (var_val = " var_val "), skipping"
continue
}
* Calculate se with explicit check
capture {
scalar se = sqrt(var_val)
if missing(se) | se <= 0 | se > 1000 {
error 498
}
}
if _rc {
di as error "Invalid se for `var' (var_val = " var_val ", se = " se ", rc = " _rc "), skipping"
continue
}
* Debug: Display scalar values
di "Estimate: " est
di "Variance: " var_val
di "Standard Error: " se
di "Values for CI calculation: est = " est ", se = " se ", invnormal(0.975) = " invnormal(0.975)
* Calculate CIs
capture scalar ci_lower = est - invnormal(0.975)*se
if _rc {
di as error "Failed to calculate ci_lower for `var' (rc = " _rc "), skipping"
continue
}
capture scalar ci_upper = est + invnormal(0.975)*se
if _rc {
di as error "Failed to calculate ci_upper for `var' (rc = " _rc "), skipping"
continue
}
* Debug: Display CI values
di "CI Lower: " ci_lower
di "CI Upper: " ci_upper
* Create matrix for collect get
matrix results = (est, ci_lower, ci_upper)
matrix colnames results = estimate ci_lower ci_upper
matrix rownames results = `coef'
* Add to collection
capture collect get results, name(`outcome'_results) tags(coef[`coef'])
if _rc {
di as error "collect get failed for `var' (rc = " _rc "), skipping"
continue
}
local ++i
di "Incremented i: `i'"
}
* Add margins results
capture noisily margins i.p203#i.period, contrast(nowald)
if _rc {
di as error "Margins failed for `outcome' (rc = " _rc "), skipping"
}
else {
* Collect margins results
capture collect get r(table), name(`outcome'_results) tags(result[margins])
if _rc {
di as error "collect get for margins failed (rc = " _rc "), skipping"
}
}
* Format collection
capture {
collect layout (coef result) (cmdset), name(`outcome'_results)
collect style cell, nformat(%8.3f)
collect style header result, title(hide)
collect label levels result estimate "Estimate" ci_lower "Lower CI" ci_upper "Upper CI" margins "Margins", modify
collect label levels coef `coef_labels', modify
collect label levels cmdset 1 "Intervention#Endline", modify
}
if _rc {
di as error "collect formatting failed (rc = " _rc ")"
}
* Preview results
di _n "Previewing collection for `outcome'"
capture collect preview, name(`outcome'_results)
if _rc {
di as error "collect preview failed (rc = " _rc ")"
}
* Save collection
capture collect export `outcome'_results, as(html) replace
if _rc {
di as error "collect export failed (rc = " _rc ")"
}
else {
di "Collection `outcome'_results saved"
}
}
Comment