*! version 1.2.8 14jun2018 Ben Jann program pshare, eclass properties(svyb svyj) version 11 if replay() { Display `0' exit } gettoken subcmd cmdline: 0, parse(", ") local length = length(`"`subcmd'"') if `"`subcmd'"'==substr("histogram",1,max(`length',1)) { capt n Histogram`cmdline' if _rc { di as err "error executing pshare histogram" exit _rc } exit } if `"`subcmd'"'==substr("stack",1,max(`length',1)) { capt n Stack`cmdline' if _rc { di as err "error executing pshare stack" exit _rc } exit } if `"`subcmd'"'==substr("contrast",1,max(`length',1)) { capt n Contrast`cmdline' // returns diopts, stack, stack2, histogram, histogram2 if _rc { di as err "error executing pshare contrast" exit _rc } ereturn local cmdline0 `"`e(cmdline)'"' ereturn local cmdline `"pshare `0'"' Display, `cdiopts' if "`stack'"!="" { Stack, `stack2' } if "`histogram'"!="" { Histogram, `histogram2' } exit } if `"`subcmd'"'!=substr("estimate",1,max(`length',1)) { capture confirm variable `subcmd' if _rc==0 { local cmdline `" `subcmd' `cmdline'"' local subcmd "estimate" local length = length(`"`subcmd'"') } } if `"`subcmd'"'==substr("estimate",1,max(`length',1)) { local version : di "version " string(_caller()) ":" Parse_opts `0' // returns lhs, options, cluster, stack, stack2, histogram, histogram2 `version' _vce_parserun pshare, mark(over) bootopts(`cluster') /// jkopts(`cluster') wtypes(pw iw fw) noeqlist: `lhs', nose `options' if "`s(exit)'" != "" { ereturn local cmdline `"pshare `0'"' if "`stack'"!="" { Stack, `stack2' } if "`histogram'"!="" { Histogram, `histogram2' } exit } Estimate`cmdline' // returns diopts, contrast, contrast2, stack, stack2, histogram, histogram2 ereturn local cmdline `"pshare `0'"' if "`contrast'"!="" { capt n Contrast_opt `contrast2' if _rc { di as err "error in contrast() option" exit _rc } } Display, `diopts' if "`stack'"!="" { Stack, `stack2' } if "`histogram'"!="" { Histogram, `histogram2' } exit } di as err `"`subcmd' invalid subcommand"' exit 198 end program Parse_opts _parse comma lhs 0 : 0 syntax [, nose CLuster(passthru) /// Stack Stack2(str asis) Histogram Histogram2(str asis) * ] c_local lhs `"`lhs'"' c_local options `options' c_local cluster `cluster' if `"`stack2'"'!="" local stack stack c_local stack `stack' c_local stack2 `"`stack2'"' if `"`histogram2'"'!="" local histogram histogram c_local histogram `histogram' c_local histogram2 `"`histogram2'"' end program Display syntax [, noHEader noTABle noGTABle Level(passthru) * ] _get_diopts diopts, `options' if `"`e(cmd)'"'!="pshare" { di as err "last pshare results not found" exit 301 } if `"`level'"'=="" local level level(`e(level)') if `"`header'"'=="" { if `"`e(contrast)'"'=="" & c(stata_version)>=12 { nobreak { Display_ereturn local cmd "total" // mimick header of -total- capture noisily break { _coef_table_header } Display_ereturn local cmd "pshare" if _rc exit _rc } } else _coef_table_header if `"`e(over)'"'!="" { _svy_summarize_legend } else di "" } if `"`table'"'=="" { if `"`e(contrast)'"'=="" & c(stata_version)>=12 { if c(stata_version)>=14 { eret di, nopvalue `level' `options' } else if c(stata_version)>=13 { quietly update if r(inst_ado)>=d(26jun2014) { eret di, nopvalue `level' `options' } else { _coef_table, cionly `level' `options' } } else { _coef_table, cionly `level' `options' } } else eret di, `level' `options' if `"`e(norm)'`e(normpop)'"'!="" { di as txt "(shares normalized with respect to total" _c if `"`e(norm)'"'!="" { di as txt " of " as res `"`e(norm)'"' _c } if `"`e(normpop)'"'!="" { if `"`e(normpop)'"'=="total" { di as txt " across all subpopulations" _c } else { di as txt `" for `e(normpop)'"' _c } } di as txt ")" } if `"`e(pvar)'"'!="" { di as txt "(percentile groups with respect to " as res `"`e(pvar)'"' as txt ")" } if `"`e(contrast)'"'!="" { di as txt "(contrasts with respect to " _c if `"`e(baseval)'"'=="+" { if `"`e(over)'"'!="" di as txt "preceding subpopulation)" else di as txt "preceding outcome variable)" } else { if `"`e(over)'"'!="" { if `"`e(baseval)'"'=="total" { di as txt `"total population)"' } else { di as txt `"`e(over)' = `e(baseval)')"' } } else di as res `"`e(baseval)'"' as txt ")" } } } if `"`gtable'"'=="" { if "`e(gini)'"!="" { matlist e(G), border(top bottom) } } end prog Display_ereturn, eclass ereturn `0' end program Histogram syntax [, /// VERTical HORizontal PROPortion NOGini Gini(str) /// keep(str) BYOPTs(str asis) LABels(str asis) BYLabels(str asis) /// Level(passthru) noci CI2(str) CIOPTs(str asis) CIBelow /// step ADDStep ADDStep2(str asis) /// OVERlay psep PSEP2(str asis) SPIKEs SPIKEs2(numlist int >0 max=1) /// prange(numlist min=2 max=2 >=0 <=100 ascending) /// base(numlist max=1) max(str asis) min(str asis) /// recast(passthru) LEGend(passthru) fxsize(passthru) fysize(passthru) /// PSTYle(passthru) PCYCle(int 15) * ] // options if `"`labels'"'!="" local bylabels `"`labels'"' if `"`addstep2'"'!="" local addstep addstep if `"`psep2'"'!="" local psep psep if "`spikes'"!="" local spikes 100 if "`spikes2'"!="" local spikes `spikes2' if "`overlay'`spikes'"!="" { local ci noci } if "`overlay'"!="" & "`psep'"!="" { di as err "overlay and psep not both allowed" exit 198 } if `"`e(cmd)'"'!="pshare" { di as err "last pshare results not found" exit 301 } if `"`recast'"'!="" { di as err "recast() not allowed" exit 198 } if "`vertical'"!="" & "`horizontal'"!="" { di as err "vertical and horizontal not both allowed" exit 198 } if `"`gini'"'=="" local gini "%9.3g" if `pcycle'<1 { di as err "pcycle() must be a positive integer" exit 198 } _get_gropts, graphopts(`options') gettwoway getallowed(plot addplot) local twowayopts `s(twowayopts)' `legend' `fxsize' `fysize' local legend local options `pstyle' `s(graphopts)' local plot `"`s(plot)'"' local addplot `"`s(addplot)'"' if `"`level'"'!="" { Parse_level, `level' // returns cilevel } else { local cilevel `e(level)' if `"`cilevel'"'=="" { // get default cilevel if e(level) not set Parse_level } } Parse_ci, `ci2' // returns citype if `"`min'"'!="" { Parse_minmax min `min' // returns min, minfmt, minnolab, minopts } if `"`max'"'!="" { Parse_minmax max `max' // returns max, maxfmt, maxnolab, maxopts } if "`min'"!="" & "`max'"!="" { if `min'>=`max' { di as err "min() must be smaller than max()" exit 198 } } // prepare CIs if "`ci'"=="" { if inlist("`citype'", "bc", "bca", "percentile") { if `"`level'"'!="" { di as err "level() not allowed with ci(`citype')" exit 198 } confirm matrix e(ci_`citype') } else if "`citype'"=="normal" { capt confirm matrix e(V) if _rc { di as err "matrix e(V) not found; cannot compute confidence intervals" exit 111 } } else { capt confirm matrix e(V) if _rc==0 { local citype normal } } } else { local citype local cibelow } // collect bars and CIs tempname S tmp0 tmp matrix `S' = e(b)' local eqs: roweq `S' local eqs: list uniq eqs foreach m in ll mid { *if "`citype'"=="" & "`m'"=="mid" continue capt mat drop `tmp' matrix `tmp0' = e(`m') foreach eq of local eqs { matrix `tmp' = nullmat(`tmp'), `tmp0'[1,"`eq':"] } matrix `S' = `S', `tmp'' } if "`citype'"!="" { if inlist("`citype'", "bc", "bca", "percentile") { local cilab " `citype'" if "`citype'"=="bc" local cilab " bias-corrected" else if "`citype'"=="bca" local cilab " BCa" matrix `S' = `S', e(ci_`citype')' } else { local cilab mata: st_matrix("`tmp'", (((1,-1)\(1,1)) * /// (st_matrix("e(b)") \ invnormal(0.5 + `cilevel'/200) * /// sqrt(diagonal(st_matrix("e(V)"))')))') mat `S' = `S', `tmp' } } // prepare keep() if `"`e(over)'"'!="" { local values `e(over_namelist)' foreach eq of local eqs { // collect over-labels from e() if "`eq'"=="total" continue local pos: list posof "`eq'" in values local overlab_`eq': word `pos' of `e(over_labels)' } } if `"`keep'"'!="" { local keeplist foreach k of local keep { if substr(`"`k'"', 1, 1)=="#" { local i = substr(`"`k'"', 2, .) capt confirm integer number `i' if _rc { di as err `"keep(): `k' invalid"' exit 198 } local i: word `i' of `eqs' } else { local i `"`k'"' if `"`e(over)'"'!="" & substr("total", 1, strlen(`"`k'"'))==`"`k'"' { local i total } if `: list i in eqs'==0 local i } if `"`i'"'=="" { di as err `"keep(): `k' invalid"' exit 198 } local keeplist `keeplist' `i' } local eqs: copy local keeplist } // collect equations, select prange, and add topcode for last bar if "`citype'"!="" local cicols ", ., ." else local cicols local bins = e(bins) if "`prange'"!="" { local prange_min: word 1 of `prange' local prange_max: word 2 of `prange' matrix `tmp' = `S'[1..`bins', 2] local bin_min 1 local bin_max 0 local pmin 0 forv i = 1/`bins' { local bin_ll = `tmp'[`i',1] if `bin_ll'<`prange_min' local bin_min = `bin_min' + 1 if `bin_min'==`i' local pmin `bin_ll' if `bin_ll'<=`prange_max' { local bin_max = `i' - 1 local pmax `bin_ll' } } if `prange_max'==100 { local bin_max `bins' local pmax 100 } local bins = `bin_max' - `bin_min' + 1 if `bins'<1 { di as err "invalid prange(): no percentile groups within specified range" exit 198 } } else { local pmin 0 local pmax 100 } matrix rename `S' `tmp0', replace tempname tmp1 mat `tmp1' = J(`bins'+1, 1, `bins'+1) forv i = 1/`bins' { mat `tmp1'[`i',1] = `i' } local i 0 foreach eq of local eqs { local ++i matrix `tmp' = `tmp0'["`eq':", 1...] if "`prange'"!="" { matrix `tmp' = `tmp'[`bin_min'..`bin_max', 1...] } matrix `tmp' = `tmp' \ (`tmp'[rowsof(`tmp'), 1], `pmax', .`cicols') matrix `tmp' = `tmp', J(rowsof(`tmp'), 1, `i'), `tmp1' matrix `S' = nullmat(`S') \ `tmp' } local keq `i' // expand data if needed and save as variables tempvar Y X MID BY PID if "`spikes'"!="" { local r = `spikes'*`keq' if _N<`r' { preserve qui set obs `r' } qui gen `Y' = . qui gen `X' = . qui gen `MID' = . qui gen `BY' = . qui gen `PID' = . local i 0 forv k = 1/`keq' { local l = (`k'-1)*(`bins'+1) + 1 local lastMID . forv j = 1/`spikes' { local ++i local p = `j' * 100/`spikes' - 50/`spikes' if `p'<`pmin' | `p'>=`pmax' continue while (`S'[`l'+1,2]<`p') { local ++l } qui replace `Y' = `S'[`l',1] in `i' qui replace `X' = `p' in `i' local newMID = `S'[`l',3] if `lastMID'!=`newMID' { qui replace `MID' = `S'[`l',3] in `i' local lastMID: copy local newMID } qui replace `BY' = `S'[`l',4] in `i' qui replace `PID' = `S'[`l',5] in `i' } } } else { local r = rowsof(`S') if _N<`r' { preserve qui set obs `r' } if "`citype'"!="" { tempvar LL UL matrix colnames `S' = `Y' `X' `MID' `LL' `UL' `BY' `PID' } else { matrix colnames `S' = `Y' `X' `MID' `BY' `PID' } svmat `S', names(col) } if "`e(pvar)'"!="" local pvar " (ordered by `e(pvar)')" if "`proportion'"!="" { qui replace `X' = `X'/100 local xti "population proportion`pvar'" qui replace `MID' = `MID'/100 } else local xti "population percentage`pvar'" lab var `X' "`xti'" if `"`e(contrast)'"'!="" { if "`e(ratio)'"!="" local yti "ratio of outcome shares" else if "`e(lnratio)'"!="" local yti "ln(ratio) of outcome shares" else local yti "difference in outcome shares" } else local yti "outcome share" lab var `Y' "`yti'" // min/max if "`min'"!="" { tempvar Ymin qui gen `Ymin' = `Y' if `Y'<`min' in 1/`r' format `Ymin' `minfmt' if "`citype'"!="" { qui replace `LL' = . if `Y'<`min' in 1/`r' qui replace `UL' = . if `Y'<`min' in 1/`r' } qui replace `Y' = `min' if `Y'<`min' in 1/`r' } if "`max'"!="" { tempvar Ymax qui gen `Ymax' = `Y' if `Y'>`max' & `Y'<. in 1/`r' format `Ymax' `maxfmt' if "`citype'"!="" { qui replace `LL' = . if `Y'>`max' in 1/`r' qui replace `UL' = . if `Y'>`max' in 1/`r' } qui replace `Y' = `max' if `Y'>`max' & `Y'<. in 1/`r' } // collect labels local hasgini = (`"`e(gini)'"'!="") & ("`nogini'"=="") if `hasgini' { tempname G matrix `G' = e(G) } local i 0 foreach eq of local eqs { local ++i gettoken bylab bylabels : bylabels if `"`bylab'"'=="" { if `"`e(over)'"'=="" { capt local bylab: var lab `eq' if `"`bylab'"'=="" { local bylab `"`eq'"' } } else { if "`eq'"=="total" local bylab "Total" else local bylab `"`overlab_`eq''"' } } if `hasgini' { local g = rownumb(`G', "`eq'") local g `:di `gini' `G'[`g',1]' local bylab `"`bylab' (Gini = `g')"' } lab def `BY' `i' `"`bylab'"', add } lab val `BY' `BY' // graph if "`overlay'"!="" { local byopt } else if `keq'>1 { local byopt if "`psep'"!="" { if "`horizontal'"!="" local byopt legend(pos(3)) else local byopt legend(pos(6)) } else if "`citype'"=="" local byopt legend(off) local byopt by(`BY', note("") `byopt' `byopts') } else { if `hasgini' { if "`eqs'"=="_" local eqs `"`e(depvar)'"' local g = rownumb(`G', "`eqs'") local g `:di `gini' `G'[`g',1]' local byopt `"subtitle("Gini = `g'")"' } else local byopt } if "`e(ratio)'"!="" & `"`base'"'=="" local base 1 if "`horizontal'"=="" { local XY `Y' `X' local cstyle c(J) if "`base'"!="" local origin yscale(range(`base')) else local origin yscale(range(0)) local base base(`base') local xtitle xtitle("`xti'") local ytitle ytitle("`yti' (`e(type)')") } else { local XY `X' `Y' local cstyle c(stepstair) if "`base'"!="" local origin xscale(range(`base')) else local origin xscale(range(0)) local xtitle xtitle("`yti' (`e(type)')") local ytitle ytitle("`xti'") } local mainplot if "`overlay'"!="" { local i 0 foreach eq of local eqs { local ++i Parse_popts o `i', `options' // returns oopts_i, options Parse_oopts `i', `oopts_`i'' // returns step_i, addstep_i addstep2_i oopts_i if "`step_`i''"=="" local step_`i' `step' if "`addstep_`i''"=="" local addstep_`i' `addstep' } local legend local i 0 local j 0 foreach eq of local eqs { local ++i local ++j local k = mod(`i'-1, `pcycle') + 1 if "`step_`i''"!="" { local mainplot `mainplot' /// (line `XY' if `BY'==`i', `cstyle' /// pstyle(p`k'line) `options' `oopts_`i'') } else if "`spikes'"!="" { local mainplot `mainplot' /// (spike `Y' `X' if `BY'==`i', `base' /// `horizontal' pstyle(p`k') `options' `oopts_`i'') } else { local mainplot `mainplot' /// (bar `Y' `X' if `BY'==`i', bartype(spanning) `base' /// `horizontal' pstyle(p`k'bar) `options' `oopts_`i'') } local legend `legend' `j' `"`:lab `BY' `i''"' if "`addstep_`i''"!="" { local ++j local mainplot `mainplot' /// (line `XY' if `BY'==`i', `cstyle' /// pstyle(p`k'line) `addstep2' `addstep2_`i'') } } local legend legend(order(`legend') all on) } else if "`psep'"!="" { // legend keys matrix `tmp' = e(ll)', e(ul)' local legend forv i = 1 / `bins' { local j `i' if "`cibelow'"!="" local j = `j' + `bins' gettoken lab psep2 : psep2 if `"`lab'"'=="" { local ll `: di %9.3g `tmp'[`i',1]' local ul `: di %9.3g `tmp'[`i',2]' local lab `ll'-`ul' } if "`horizontal'"!="" { local legend `j' `"`lab'"' `legend' } else { local legend `legend' `j' `"`lab'"' } } if "`horizontal'"!="" { local legend legend(order(`legend') all stack keygap(*.5) rowgap(*1.5) cols(1) pos(3)) } else { local legend legend(order(`legend') all stack rows(1) pos(6)) } // p#() forv i = 1/`bins' { Parse_popts p `i', `options' // returns popts_i, options Parse_ciopts p `i', `popts_`i'' // returns pciopts_i, ppsty_i, popts_i } // plots local i 0 forv i = 1/`bins' { local k = mod(`i'-1, `pcycle') + 1 if "`step'"!="" { local mainplot `mainplot' /// (line `XY' if inlist(`PID',`i',`i'+1), `cstyle' /// pstyle(p`k'line) `options' `ppsty_`i'' `popts_`i'') } else if "`spikes'"!="" { local mainplot `mainplot' /// (spike `Y' `X' if `PID'==`i', `base' /// `horizontal' pstyle(p`k') `options' `ppsty_`i'' `popts_`i'') } else { local mainplot `mainplot' /// (bar `Y' `X' if inlist(`PID',`i',`i'+1), bartype(spanning) `base' /// `horizontal' pstyle(p`k'bar) `options' `ppsty_`i'' `popts_`i'') } if "`citype'"!="" { local ciplot `ciplot' (rcap `LL' `UL' `MID' /// if `PID'==`i', `horizontal' /// pstyle(p`k') `pstyle' `ciopts' `ppsty_`i'' `pciopts_`i'') } } if "`addstep'"!="" { local mainplot `mainplot' (line `XY', `cstyle' /// pstyle(p1line) `addstep2') } } else { local j 2 if "`step'"!="" { local mainplot (line `XY', `cstyle' /// pstyle(p1line) `options') } else if "`spikes'"!="" { local mainplot (spike `Y' `X', `base' /// `horizontal' pstyle(p1) `options' ) } else { local mainplot (bar `Y' `X', bartype(spanning) `base' /// `horizontal' pstyle(p1bar) `options' ) } if "`addstep'"!="" { local ++j local mainplot `mainplot' (line `XY', `cstyle' /// pstyle(p1line) `addstep2') } if "`citype'"!="" { if "`cibelow'"!="" /// local legend legend(order(`j' 1 "`cilevel'%`cilab' CI")) else /// local legend legend(order(1 `j' "`cilevel'%`cilab' CI")) local ciplot (rcap `LL' `UL' `MID', `horizontal' pstyle(p2) `ciopts') } else if `"`byopt'"'=="" local legend legend(off) } if "`citype'"!="" { if "`cibelow'"!="" { local mainplot `ciplot' `mainplot' } else { local mainplot `mainplot' `ciplot' } } if "`min'"!="" & "`minnolab'"=="" { if "`horizontal'"=="" { local mainplot `mainplot' (scatter `Y' `MID' if `Ymin'<., ms(i) /// mlabel(`Ymin') mlabc(fg) mlabp(6) `minopts') } else { local mainplot `mainplot' (scatter `MID' `Y' if `Ymin'<., ms(i) /// mlabel(`Ymin') mlabc(fg) mlabp(9) `minopts') } } if "`max'"!="" & "`maxnolab'"=="" { if "`horizontal'"=="" { local mainplot `mainplot' (scatter `Y' `MID' if `Ymax'<., ms(i) /// mlabel(`Ymax') mlabc(fg) mlabp(12) `maxopts') } else { local mainplot `mainplot' (scatter `MID' `Y' if `Ymax'<., ms(i) /// mlabel(`Ymax') mlabc(fg) mlabp(3) `maxopts') } } twoway `mainplot' || `plot' || `addplot' || in 1/`r', /// `xtitle' `ytitle' `origin' `legend' `twowayopts' `byopt' end program Parse_level syntax [, Level(cilevel) ] c_local cilevel `"`level'"' end program Parse_ci syntax [, NORmal bc bca Percentile ] local citype `normal' `bc' `bca' `percentile' if `:list sizeof citype'>1 { di as err "only one of normal, bc, bca, and percentile allowed in ci()" exit 198 } c_local citype `citype' end program Parse_minmax gettoken minmax 0 : 0 syntax anything(id="`minmax'(): number") [, Format(str) NOLABels * ] confirm number `anything' if `"`format'"'=="" local format "%9.3g" confirm numeric format `format' c_local `minmax' `anything' c_local `minmax'fmt `format' c_local `minmax'nolab `nolabels' c_local `minmax'opts `options' end program Stack syntax [ , /// VERTical HORizontal REVerse /// keep(str) Sort Sort2(str) /// Values Values2(str) NOGini Gini(str) /// LABels(str asis) PLABels(str asis) /// BARWidth(real 0.75) /// LEGend(passthru) fxsize(passthru) fysize(passthru) PCYCle(int 15) * /// ] // options if `"`e(cmd)'"'!="pshare" { di as err "last pshare results not found" exit 301 } if "`vertical'"!="" & "`horizontal'"!="" { di as err "vertical and horizontal not both allowed" exit 198 } if `pcycle'<1 { di as err "pcycle() must be a positive integer" exit 198 } if "`vertical'`horizontal'"=="" local horizontal horizontal if `"`values2'"'!="" local values values else if "`values'"!="" local values2 "%9.3g" if "`values'"!="" local valuespos mlabposition(0) _get_gropts, graphopts(`options') gettwoway getallowed(plot addplot) local twowayopts `s(twowayopts)' `legend' `fxsize' `fysize' local legend local options `s(graphopts)' local plot `"`s(plot)'"' local addplot `"`s(addplot)'"' if `"`sort'`sort2'"'!="" { Parse_sort, `sort2' // returns sort, sortgini, descending, tlast, tfirst if `"`over'"'!="" & "`tlast'`tfirst'"!="" local total total } if `"`gini'"'=="" local gini "%9.3g" // collect results local bins = e(bins) local k_eq = e(k_eq) tempname S matrix `S' = e(b)' local eqs: roweq `S' local eqs: list uniq eqs local hasgini = (`"`e(gini)'"'!="") if `hasgini' { tempname G matrix `G' = e(G) } else if "`sortgini'"!="" { di as err "Gini coefficients not found; apply the gini option to pshare estimate" exit 111 } // keep() local eqs0: copy local eqs if `"`keep'"'!="" { local keeplist foreach k of local keep { if substr(`"`k'"', 1, 1)=="#" { local i = substr(`"`k'"', 2, .) capt confirm integer number `i' if _rc { di as err `"keep(): `k' invalid"' exit 198 } local i: word `i' of `eqs' } else { local i `"`k'"' if `"`e(over)'"'!="" & /// substr("total", 1, strlen(`"`k'"'))==`"`k'"' { local i total } if `: list i in eqs'==0 local i } if `"`i'"'=="" { di as err `"keep(): `k' invalid"' exit 198 } local keeplist `keeplist' `i' } local eqs: copy local keeplist } // sort() if "`sort'"!="" { mata: pshare_stack_sort() // returns local eqs0 local eqs: list eqs0 & eqs } // select results and define plot positions tempname id matrix `id' = J(`bins', 1, .) forv i = 1/`bins' { matrix `id'[`i',1] = `i' } tempname S0 tmp matrix rename `S' `S0' if `hasgini' { tempname G0 matrix rename `G' `G0' } local i 0 foreach eq of local eqs { local ++i matrix `tmp' = `S0'["`eq':", 1...] matrix `tmp' = `tmp', J(rowsof(`tmp'), 1, `i'), `id' matrix `S' = nullmat(`S') \ `tmp' if `hasgini' { if "`eq'"=="_" { // happens if there is just one set matrix `G' = `G0' } else { matrix `G' = nullmat(`G') \ `G0'["`eq'", 1...] } } } local nstacks `i' matrix `S' = J(rowsof(`S'), 2, .), `S' local reverse = "`reverse'"!="" if `reverse' local i = `nstacks'*`bins'+1 else local i 0 forv j = 1/`nstacks' { local neg0 0 local pos0 0 forv k = 1/`bins' { local i = `i' + cond(`reverse', -1, 1) if (`S'[`i', 3])<0 local dir neg else local dir pos matrix `S'[`i', 1] = ``dir'0' // LL local `dir'0 = ``dir'0' + `S'[`i', 3] matrix `S'[`i', 2] = ``dir'0' // UL } } // expand data if needed and save results to data local r = rowsof(`S') if _N<`r' { preserve qui set obs `r' } tempvar LL UL Y X ID MID matrix colnames `S' = `LL' `UL' `Y' `X' `ID' svmat `S', names(col) qui gen `MID' = `LL' + (`UL' - `LL')/2 if "`values'"!="" { format `Y' `values2' } // labels local ylab local i 0 local over_namelist `e(over_namelist)' foreach eq of local eqs { local ++i gettoken lab labels : labels if `"`lab'"'=="" { if `"`e(over)'"'=="" { local eqi "`eq'" if "`eqi'"=="_" local eqi "`e(depvar)'" capt local lab: var lab `eqi' if `"`lab'"'=="" { local lab "`eqi'" } } else { if "`eq'"=="total" local lab "Total" else { local pos: list posof "`eq'" in over_namelist local lab: word `pos' of `e(over_labels)' } } } lab def `X' `i' `"`lab'"', add } lab val `X' `X' // legend keys matrix `tmp' = e(ll)', e(ul)' local legend forv i = 1 / `bins' { gettoken lab plabels : plabels if `"`lab'"'=="" { if `reverse' { local j = `bins' + 1 - `i' local ll `: di %9.3g 100-`tmp'[`j',2]' local ul `: di %9.3g 100-`tmp'[`j',1]' } else { local ll `: di %9.3g `tmp'[`i',1]' local ul `: di %9.3g `tmp'[`i',2]' } local lab `ll'-`ul' } if "`horizontal'"!="" { local legend `legend' `i' `"`lab'"' } else { local legend `i' `"`lab'"' `legend' } } // gini labels if `hasgini' & "`nogini'"=="" { local ginilbls forv i = 1 / `nstacks' { local lab `: di `gini' `G'[`i',1]' local ginilbls `ginilbls' `i' "`lab'" } if "`horizontal'"!="" { local giniaxis yaxis(1 2) local ginilbls ylabels(`ginilbls', axis(2) notick labgap(*2) /// labsize(*.8) angle(horizontal)) ytitle("Gini", axis(2) /// size(*.8) orient(rvertical)) } else { local giniaxis xaxis(1 2) local ginilbls xlabels(`ginilbls', axis(2) notick labgap(*2) /// labsize(*.8)) xtitle("Gini", axis(2) size(*.8)) } } // plots local plots local vplots forv i = 1 / `bins' { Parse_popts p `i', `options' // returns popts_i, options } forv i = 1 / `bins' { if `reverse' local j = `bins' + 1 - `i' else local j `i' local plots `plots' (rbar `LL' `UL' `X' if `ID'==`j', /// barwidth(`barwidth') `giniaxis' `options' `popts_`i'' `vertical'`horizontal') if "`values'"!="" { local k = mod(`i'-1, `pcycle') + 1 if "`horizontal'"!="" { local vplots `vplots' (scatter `X' `MID' if `ID'==`j', /// mlabel(`Y') ms(i) psty(p`k') `valuespos' `options' `popts_`i'') } else { local vplots `vplots' (scatter `MID' `X' if `ID'==`j', /// mlabel(`Y') ms(i) psty(p`k') `valuespos' `options' `popts_`i'') } } } // draw graph local xti `"`e(type)'"' if "`horizontal'"!="" { local yscale = `nstacks'+0.5 local yscale yscale(range(0.5 `yscale') reverse) two `plots' `vplots' || `plot' || `addplot' /// ||, plotregion(margin(zero)) /// ylabels(1(1)`nstacks', val notick angle(horizontal) labgap(*2) nogrid) /// yti("") `yscale' /// xlabel(, grid) xtitle("outcome share (`xti')") /// legend(order(`legend') all stack rows(1) pos(6)) /// `ginilbls' pcycle(`pcycle') `twowayopts' } else { local xscale = `nstacks'+0.5 local xscale xscale(range(0.5 `xscale')) two `plots' `vplots' || `plot' || `addplot' /// ||, plotregion(margin(zero)) /// xlabels(1(1)`nstacks', val notick labgap(*2) nogrid) xti("") `xscale' /// ytitle("outcome share (`xti')") ylabel(, grid angle(horizontal)) /// legend(order(`legend') all stack keygap(*.5) rowgap(*1.5) cols(1) pos(3)) /// `ginilbls' `twowayopts' } end program Parse_sort syntax [, Gini Descending TLast TFirst ] if "`tlast'"!="" & "`tfirst'"!="" { di as err "sort(): only one of tfirst and tlast allowed" exit 198 } c_local sort sort c_local sortgini `gini' c_local descending `descending' c_local tlast `tlast' c_local tfirst `tfirst' end program Parse_popts gettoken p 0 : 0, parse(" ,") gettoken j 0 : 0, parse(" ,") syntax [, `p'`j'(str asis) * ] c_local `p'opts_`j' ``p'`j'' c_local options `options' end program Parse_oopts gettoken j 0 : 0, parse(" ,") syntax [, step ADDStep ADDStep2(str asis) * ] if `"`addstep2'"'!="" local addstep addstep c_local step_`j' `step' c_local addstep_`j' `addstep' c_local addstep2_`j' `"`addstep2'"' c_local oopts_`j' `options' end program Parse_ciopts gettoken p 0 : 0, parse(" ,") gettoken j 0 : 0, parse(" ,") syntax [, CIOPTs(str asis) PSTYle(passthru) * ] c_local `p'ciopts_`j' `ciopts' c_local `p'psty_`j' `pstyle' c_local `p'opts_`j' `options' end program Estimate, sortpreserve eclass // syntax syntax varlist(numeric) [if] [in] [pw iw fw/], [ /// PRoportion percent DENsity sum AVErage GENERALized /// NORMalize(str) base(str) GRANDtotal gini /// Nquantiles(numlist int max=1 >=2) /// Percentiles(numlist >=0 <=100 ascending) /// pvar(varname numeric) step /// over(varname numeric) Total /// Contrast Contrast2(str) /// Stack Stack2(str asis) Histogram Histogram2(str asis) /// vce(str) CLuster(varname) svy SVY2(str) nose /// /// display options Level(cilevel) noHEader noTABle noGTABle * /// ] // varlist local varlist: list uniq varlist // remove repeated varnames local ndepv: list sizeof varlist // type of results local rtype `proportion' `percent' `density' `sum' `average' `generalized' if `: list sizeof rtype'>1 { local rtype: list retok rtype di as err "`rtype': only one allowed" exit 198 } if "`rtype'"=="" local rtype proportion // normalize() if `"`base'"'!="" & `"`normalize'"'=="" { // backward compatability local normalize: copy local base } if `"`normalize'"'!="" { Parse_normalize, norm(`normalize') // returns norm, normnum, normpop if "`norm'"=="*" local norm `varlist' } if "`grandtotal'"!="" & `"`normpop'"'=="" local normpop total // backward compatability if inlist("`rtype'", "sum", "average", "generalized") { if `"`norm'`normnum'`normpop'"'!="" { di as err "normalize() not allowed with `rtype'" exit 198 } } // percentile cutoffs, labels for percentile bins if "`nquantiles'"!="" & "`percentiles'"!="" { di as err "nquantiles() and percentiles() not both allowed" exit 198 } if "`nquantiles'"=="" local nquantiles 5 if "`percentiles'"!="" { local tmp foreach p of local percentiles { // get rid of 0 and 100 if specified if inlist(`p', 0, 100) continue local tmp `tmp' `p' } local percentiles `tmp' if "`percentiles'"=="" { di as err "percentiles() invalid -- too few elements specified" exit 122 } } else { forv i = 1/`=`nquantiles'-1' { local percentiles `percentiles' `=`i'*100/`nquantiles'' } } local bins 0 local p_ll 0 local binlbls foreach p_ul in `percentiles' 100 { local ++bins local ll `: di %9.5g `p_ll'' local ul `: di %9.5g `p_ul'' local binlbls `binlbls' `ll'-`ul' local p_ll `p_ul' } local binlbls: subinstr local binlbls "." ",", all // to fix problem with decimal point // over(), total if "`over'"=="" { if "`total'"!="" { di as err "total only allowed if over() is specified" exit 198 } if "`grandtotal'"!="" { di as err "grandtotal only allowed if over() is specified" exit 198 } } else if `ndepv'>1 { di as err "over() not allowed if multiple variables are specified" exit 198 } // vce local vce0 `"`vce'"' if `"`cluster'"'!="" { if "`se'"!="" { di as err "cluster() and nose not both allowed" exit 198 } if `"`svy'`svy2'"'!="" { di as err "cluster() and svy() not both allowed" exit 198 } if `"`vce'"'!="" { di as err "cluster() and vce() not both allowed" exit 198 } local vce `"cluster `cluster'"' local cluster } if `"`vce'"'!="" { if "`se'"!="" { di as err "vce() and nose not both allowed" exit 198 } if `"`svy'`svy2'"'!="" { di as err "vce() and svy() not both allowed" exit 198 } gettoken vce clustvar : vce if `"`vce'"'=="analytic" & `"`clustvar'"'=="" { local vce analytic } else if `:list sizeof clustvar'==1 & /// substr("cluster", 1, max(2, strlen(`"`vce'"')))==`"`vce'"' { local vce cluster local vcetype Robust gettoken clustvar : clustvar local clustopt cluster(`clustvar') } else { di as err "invalid vce()" exit 198 } } else if "`se'"=="" { local vce analytic } // svy if `"`svy'`svy2'"'!="" { if "`weight'"!="" { di as err "weights not allowed with svy; supply weights to {help svyset}" exit 101 } if "`se'"!="" { di as err "svy() and nose not both allowed" exit 198 } local svy svy if `"`svy2'"'!="" { local svy2 `"subpop(`svy2')"' } } // display options local levelopt level(`level') _get_diopts diopts, `options' c_local diopts `header' `table' `gtable' `levelopt' `diopts' // graph options if `"`stack2'"'!="" local stack stack c_local stack `stack' c_local stack2 `"`stack2'"' if `"`histogram2'"'!="" local histogram histogram c_local histogram `histogram' c_local histogram2 `"`histogram2'"' // sample and weights marksample touse markout `touse' `clustvar' `over' `pvar' `norm' if "`svy'"!="" { if c(stata_version)>=14 { tempvar subpop exp _svy_setup `touse' `subpop' `exp', svy `svy2' local weight `"`r(wtype)'"' if "`weight'"=="" { drop `exp' local exp } } else { tempvar subpop _svy_setup `touse' `subpop', svy `svy2' local weight `"`r(wtype)'"' local exp `"`r(wvar)'"' } if `"`r(vce)'"'!="linearized" { di as err "option svy is only allowed if VCE is set to linearized; " /// `"use the {helpb svy} prefix command for `r(vce)' survey estimation"' exit 498 } local svy_posts `"`r(poststrata)'"' local svy_postw `"`r(postweight)'"' if `"`svy_posts'"'!="" { if "`weight'"!="" local wexp `"[`weight' = `exp']"' tempvar exp svygen post double `exp' `wexp' if `touse', /// posts(`svy_posts') postw(`svy_postw') if "`weight'"=="" local weight pweight } } else local subpop `touse' if "`over'"!="" { capt assert ((`over'==floor(`over')) & (`over'>=0)) if `subpop' if _rc { di as err "variable in over() must be integer and nonnegative" exit 452 } qui levelsof `over' if `subpop', local(overvals) local N_over: list sizeof overvals local over_labels foreach overval of local overvals { local over_labels `"`over_labels' `"`: label (`over') `overval''"'"' } local over_labels: list clean over_labels if `"`normpop'"'!="" & `"`normpop'"'!="total" { if substr(`"`normpop'"', 1, 1)=="#" { local normpop = substr(`"`normpop'"', 2, .) capt confirm integer number `normpop' if _rc { di as err `"normalize(): #`normpop' invalid"' exit 198 } local normpop: word `normpop' of `overvals' } else { capt confirm integer number `normpop' if _rc { di as err `"normalize(): `normpop' invalid"' exit 198 } local normpop: list overvals & normpop } if `"`normpop'"'=="" { di as err "normalize(): specified subpopulation not found" exit 198 } } } else if `"`normpop'"'!="" { di as err "reference population in normalize() only allowed if over() is specified" exit 198 } if "`weight'"!="" { capt confirm variable `exp' if _rc { local wexp `"= `exp'"' tempvar exp qui gen double `exp' `wexp' } local wgt `"[`weight' = `exp']"' if "`weight'"=="pweight" { local swgt `"[aw = `exp']"' } else { local swgt `"`wgt'"' } } else { local wgt local swgt } if "`svy'"!="" & `"`exp'"'!="" { su `touse' if `touse', meanonly } else { su `touse' `swgt' if `touse', meanonly } local N = r(N) if `N'==0 error 2000 // contrast if `"`contrast2'"'!="" local contrast contrast if "`contrast'"!="" { if "`over'"=="" { if `ndepv'==1 { di as err "contrast not allowed with single outcome variable" exit 198 } } else { if `N_over'==1 { di as err "contrast not allowed with single subpopulation" exit 198 } } } c_local contrast `contrast' c_local contrast2 `"`contrast2'"' // check whether variance estimation is supported if "`se'"=="" { if inlist("`weight'", "fweight", "iweight") { if "`vce'"=="analytic" & `"`vce0'"'=="" { di as txt "(variance estimation not supported with `weight's)" local se nose } else if "`svy'`vce'"!="" { di as err "variance estimation not supported with `weight's" exit 498 } } } // generate helper variable for uvars if "`se'"=="" & ("`rtype'"=="sum" | "`normnum'"!="") { tempvar avar if "`svy'"!="" { qui gen byte `avar' = 0 if `touse' } else { Generate_avar `avar', touse(`touse') c(`clustvar') t(`normnum') /// n(`N') w(`exp') `svy' } } // compute shares tempname S local uvars if "`pvar'"!="" tempvar vtmp if "`norm'"!="" { tempname normvar qui gen double `normvar' = 0 if `touse' foreach v of local norm { qui replace `normvar' = `normvar' + `v' if `touse' } } if "`touse'"=="`subpop'" local offset 0 else { qui count if `touse' & `subpop'==0 local offset = r(N) } if "`over'"=="" { if "`gini'"!="" { tempname G mat `G' = J(`ndepv', 1, .) mat rown `G' = `varlist' mat coln `G' = "Gini" } local k_eq 0 foreach v of local varlist { local uvarsi if "`se'"=="" { forv i = 1/`bins' { tempvar u qui gen double `u' = . local uvarsi `uvarsi' `u' } local uvars `uvars' `uvarsi' } local ++k_eq if `ndepv'>1 local eq "`v'" else local eq if "`pvar'"!="" { sort `subpop' `pvar' `v' `exp' if `"`exp'"'!="" { qui by `subpop' `pvar': gen `vtmp' = /// sum(`exp'*`v')/sum(`exp') if `subpop' } else { qui by `subpop' `pvar': gen `vtmp' = sum(`v')/_N if `subpop' } qui by `subpop' `pvar': replace `vtmp' = `vtmp'[_N] if `subpop' } else { sort `subpop' `v' `exp' local vtmp `v' } _Estimate `vtmp' `v' `wgt', `gini' binlbls(`binlbls') /// touse(`touse') subpop(`subpop') offset(`offset') /// pvar(`pvar') percentiles(`percentiles') `step' /// rtype(`rtype') eqlab(`eq') uvars(`uvarsi') /// normnum(`normnum') avar(`avar') normvar(`normvar') mat `S' = nullmat(`S') \ r(S) if "`gini'"!="" { mat `G'[`k_eq',1] = r(gini) } if "`pvar'"!="" drop `vtmp' } } else { tempname _N mat `_N' = J(`N_over'+("`total'"!=""), 1, .) mat rown `_N' = `overvals' `total' mat coln `_N' = "N" if "`gini'"!="" { tempname G mat `G' = `_N' mat coln `G' = "Gini" } local ytotal local wtotal if "`normpop'"!="" { if "`normpop'"=="total" local ifnormpop else local ifnormpop `" & `over'==`normpop'"' if "`normnum'"!="" { su `varlist' `swgt' if `subpop'`ifnormpop', meanonly tempname wtotal scalar `wtotal' = r(sum_w) } else if "`normvar'"!="" { su `normvar' `swgt' if `subpop'`ifnormpop', meanonly tempname ytotal wtotal scalar `ytotal' = r(sum) scalar `wtotal' = r(sum_w) } else { su `varlist' `swgt' if `subpop'`ifnormpop', meanonly tempname ytotal wtotal scalar `ytotal' = r(sum) scalar `wtotal' = r(sum_w) } } if "`pvar'"!="" { sort `subpop' `over' `pvar' `varlist' `exp' if `"`exp'"'!="" { qui by `subpop' `over' `pvar': gen `vtmp' = /// sum(`exp'*`varlist')/sum(`exp') if `subpop' } else { qui by `subpop' `over' `pvar': gen `vtmp' = sum(`varlist')/_N if `subpop' } qui by `subpop' `over' `pvar': replace `vtmp' = `vtmp'[_N] if `subpop' } else { sort `subpop' `over' `varlist' `exp' local vtmp `varlist' } local k_eq 0 local offseti `offset' foreach overval of local overvals { local ++k_eq local uvarsi if "`se'"=="" { forv i = 1/`bins' { tempvar u qui gen double `u' = . local uvarsi `uvarsi' `u' } local uvars `uvars' `uvarsi' } _Estimate `vtmp' `varlist' `wgt', `gini' binlbls(`binlbls') /// touse(`touse') subpop(`subpop') offset(`offseti') /// pvar(`pvar') percentiles(`percentiles') `step' /// rtype(`rtype') eqlab(`overval') uvars(`uvarsi') /// normnum(`normnum') avar(`avar') normvar(`normvar') /// normpop(`normpop') over(`over') overval(`overval') /// ytotal(`ytotal') wtotal(`wtotal') mat `_N'[`k_eq', 1] = r(N) mat `S' = nullmat(`S') \ r(S) if "`gini'"!="" { mat `G'[`k_eq',1] = r(gini) } } if "`total'"!="" { local ++k_eq if "`pvar'"!="" { drop `vtmp' sort `subpop' `pvar' `varlist' `exp' if `"`exp'"'!="" { qui by `subpop' `pvar': gen `vtmp' = /// sum(`exp'*`varlist')/sum(`exp') if `subpop' } else { qui by `subpop' `pvar': gen `vtmp' = sum(`varlist')/_N if `subpop' } qui by `subpop' `pvar': replace `vtmp' = `vtmp'[_N] if `subpop' } else { sort `subpop' `varlist' `exp' local vtmp `varlist' } local uvarsi if "`se'"=="" { forv i = 1/`bins' { tempvar u qui gen double `u' = . local uvarsi `uvarsi' `u' } local uvars `uvars' `uvarsi' } _Estimate `vtmp' `varlist' `wgt', `gini' binlbls(`binlbls') /// touse(`touse') subpop(`subpop') offset(`offset') /// pvar(`pvar') percentiles(`percentiles') `step' /// rtype(`rtype') eqlab(total) uvars(`uvarsi') /// normnum(`normnum') avar(`avar') normvar(`normvar') /// normpop(`normpop') over(`over') /// ytotal(`ytotal') wtotal(`wtotal') mat `_N'[`k_eq', 1] = r(N) mat `S' = nullmat(`S') \ r(S) if "`gini'"!="" { mat `G'[`k_eq',1] = r(gini) } } } if "`se'"=="" { tempname V if "`svy'"!="" { qui svy, `svy2': total `uvars' if `touse' } else { qui total `uvars' `wgt' if `touse', `clustopt' } matrix `V' = e(V) local rank = e(rank) local df_r = e(df_r) if "`vce'"=="cluster" { local N_clust = e(N_clust) } if "`svy'"!="" { local svy_scalars foreach l in N_sub N_strata N_strata_omit singleton census N_pop /// N_subpop N_psu N_poststrata stages { if e(`l')<. { local svy_`l' = e(`l') local svy_scalars `svy_scalars' `l' } } local svy_macros foreach l in prefix wtype wvar wexp singleunit strata psu fpc /// poststrata postweight vce vcetype mse subpop adjust { local svy_`l' `"`e(`l')'"' local svy_macros `svy_macros' `l' } forv l=1/`svy_stages' { local svy_su`l' `"`e(su`l')'"' local svy_fpc`l' `"`e(fpc`l')'"' local svy_weight`l' `"`e(weight`l')'"' local svy_strata`l' `"`e(strata`l')'"' local svy_macros `svy_macros' su`l' fpc`l' weight`l' strata`l' } local svy_matrices foreach l in V_srs V_srssub V_srswr V_srssubwr _N_strata_single /// _N_strata_certain _N_strata _N_postsum _N_postsize { capt confirm matrix e(`l') if _rc==0 { tempname svy_`l' mat `svy_`l'' = e(`l') local svy_matrices `svy_matrices' `l' } } } } tempname b mat `b' = `S'[1..., 7]' // post results if "`se'"=="" { local coln: colfullnames `b' mat coln `V' = `coln' mat rown `V' = `coln' } else if "`weight'"=="iweight" { // add empty e(V) for svy tempname V mat `V' = `b'' * `b' * 0 } eret post `b' `V' `wgt', obs(`N') esample(`touse') eret local cmd "pshare" eret local depvar "`varlist'" eret local pvar "`pvar'" eret local percentiles "`percentiles'" eret local step "`step'" eret local type "`rtype'" eret local norm "`norm'`normnum'" if "`normpop'"!="" { if "`normpop'"=="total" eret local normpop total else eret local normpop "`over' = `normpop'" } local title "Percentile shares (`rtype')" if "`svy'"!="" local title "Survey: `title'" eret local title "`title'" eret scalar k_eq = `k_eq' eret scalar bins = `bins' local i 0 tempname tmp foreach m in ll ul mid prop L_ll L_ul { local ++i mat `tmp' = `S'[1..., `i']' mat rown `tmp' = "`m'" eret mat `m' = `tmp' } if "`over'"!="" { eret local total "`total'" eret local over_labels `"`over_labels'"' eret local over_namelist `"`overvals'"' eret local over "`over'" eret scalar N_over = `N_over' eret mat _N = `_N' } if "`gini'"!="" { eret local gini "gini" eret mat G = `G' } if "`se'"=="" { eret local vcetype "`vcetype'" eret local vce "`vce'" eret scalar rank = `rank' eret scalar df_r = `df_r' eret scalar level = `level' if "`vce'"=="cluster" { eret local clustvar "`clustvar'" eret scalar N_clust = `N_clust' } if "`svy'"!="" { foreach l of local svy_scalars { eret scalar `l' = `svy_`l'' } foreach l of local svy_macros { eret local `l' `"`svy_`l''"' } foreach l of local svy_matrices { if substr("`l'", 1, 1)=="V" { mat coln `svy_`l'' = `coln' mat rown `svy_`l'' = `coln' } eret matrix `l' = `svy_`l'' } } } Return_clear // for some reason -return clear- does not delete r(S) end program Parse_normalize // get over from norm(over:total) syntax, norm(str) gettoken normpop norm : norm, parse(":") if `"`normpop'"'==":" { local normpop local norm `":`norm'"' } else if `"`normpop'"'=="." local normpop else if `"`norm'"'!="" { if substr("total", 1, strlen(`"`normpop'"'))==`"`normpop'"' { local normpop total } } if `"`norm'"'==":" { c_local norm c_local normpop `"`normpop'"' c_local normnum exit } if `"`norm'"'=="" { local norm `":`normpop'"' local normpop } // get total from norm(over:total) gettoken norm 0 : norm, parse(":") local 0 `", norm(`0')"' capt syntax, norm(numlist max=1) // norm(#) if _rc==0 { if `norm'==0 { di as err "normalize() cannot be 0" exit 198 } c_local norm c_local normpop `"`normpop'"' c_local normnum `norm' exit } syntax [, norm(str) ] if `"`norm'"'=="." local norm if `"`norm'"'=="*" | `"`norm'"'=="" { c_local norm `norm' c_local normpop `"`normpop'"' c_local normnum exit } syntax, norm(varlist numeric) c_local norm `norm' c_local normpop `"`normpop'"' c_local normnum end program Generate_avar syntax name(name=avar) [, touse(str) c(str) t(str) n(str) w(str) /*svy*/ ] if "`t'"=="" local t 1 if "`w'"=="" local w `touse' if "`c'"=="" { tempvar c qui gen double `c' = _n if `touse' } sort `touse' `c' // count number of clusters qui by `touse': gen double `avar' = sum(`c'!=`c'[_n-1]) if `touse' qui by `touse': replace `avar' = `avar'[_N] if `touse' // aggregate weights within clusters qui by `touse' `c': replace `avar' = `avar' * sum(`w') if `touse' qui by `touse' `c': replace `avar' = `avar'[_N] if `touse' // invert qui replace `avar' = `t' / `avar' if `touse' end program _Estimate, rclass syntax varlist [pw iw fw/] [, /// gini binlbls(str) touse(str) subpop(str) offset(str) /// pvar(str) percentiles(str) step rtype(str) eqlab(str) uvars(str) /// normnum(str) avar(str) normvar(str) normpop(str) /// over(str) overval(str) ytotal(str) wtotal(str) ] gettoken y y0 : varlist // in case of pvar(): y contains averaged values gettoken y0 : y0 // y0 contrains original values if ("`y'"=="`y0'") local y0 local touse0 `touse' if "`over'"=="" local touse `subpop' else { if "`overval'"=="" local touse `subpop' // total else { tempvar touse qui gen byte `touse' = `subpop' & (`over'==`overval') } if "`normpop'"!="" { if "`normpop'"=="total" local touseR `subpop' else { tempvar touseR qui gen byte `touseR' = `subpop' & (`over'==`normpop') } } } tempname S mata: pshare_compute() // sets matrix S and local N (and G), replaces uvars mat rown `S' = `binlbls' if "`eqlab'"!="" { mat roweq `S' = "`eqlab'" } return scalar N = `N' return matrix S = `S' if "`gini'"!="" { return scalar gini = `G' } c_local offseti `offset' end program Return_clear, rclass local x end program Contrast_opt syntax [anything] [, Ratio LNRatio ] Contrast `anything' , `ratio' `lnratio' end program Contrast, eclass // syntax syntax [anything(name=base)] [, Ratio LNRatio /// Stack Stack2(str asis) Histogram Histogram2(str asis) * ] if `"`e(cmd)'"'!="pshare" { di as err "last pshare results not found" exit 301 } if `"`e(contrast)'"'!="" { di as err "contrast already applied" exit 321 } if "`lnratio'"!="" & "`ratio'"!="" { di as err "only one of ratio and lnratio allowed" exit 198 } c_local cdiopts `options' if `"`stack2'"'!="" local stack stack c_local stack `stack' c_local stack2 `"`stack2'"' if `"`histogram2'"'!="" local histogram histogram c_local histogram `histogram' c_local histogram2 `"`histogram2'"' // parse base local varlist `"`e(depvar)'"' local ndepv: list sizeof varlist local over `"`e(over)'"' local N_over = e(N_over) local overvals `"`e(over_namelist)'"' local total `"`e(total)'"' if `"`over'"'=="" { if `ndepv'==1 { di as err "contrast not allowed with single outcome variable" exit 198 } } else { if `N_over'==1 { di as err "contrast not allowed with single subpopulation" exit 198 } } local baseval if `"`base'"'!="" { if `: list sizeof base'>1 { di as err `"contrast: `base' invalid"' exit 198 } if substr(`"`base'"', 1, 1)=="#" { local baseval = substr(`"`base'"', 2, .) capt confirm integer number `baseval' if _rc { di as err `"contrast: `base' invalid"' exit 198 } if `"`over'"'=="" { local baseval: word `baseval' of `varlist' } else { local baseval: word `baseval' of `overvals' } } else { if `"`over'"'=="" { capt unab baseval: `base' if _rc local baseval `"`base'"' if `: list baseval in varlist'==0 local baseval } else { local baseval `"`base'"' if `: list baseval in overvals'==0 local baseval } } if "`baseval'"=="" { di as err `"contrast: `base' invalid"' exit 498 } } if "`baseval'"=="" { if `"`total'"'!="" local baseval total else local baseval + } // transform results tempname b mat `b' = e(b) capt confirm matrix e(V) if _rc local se nose else { local se tempname V mat `V' = e(V) } local k_eq = e(k_eq) local --k_eq local bins = e(bins) if "`se'"=="" & "`ratio'`lnratio'"!="" { Ratio_V `ratio'`lnratio' `b' `V' "`baseval'" } else { local eqs: coleq `b' local eqs: list uniq eqs tempname B tmp tmp2 matrix rename `b' `tmp' if "`baseval'"!="+" { matrix `B' = `tmp'[1...,"`baseval':"] } local eq0 foreach eq of local eqs { if "`baseval'"=="`eq'" continue if "`baseval'"=="+" { if "`eq0'"=="" { local eq0 "`eq'" continue } else { matrix `B' = `tmp'[1...,"`eq0':"] local eq0 "`eq'" } } matrix `tmp2' = `tmp'[1...,"`eq':"] if "`ratio'"!="" { mata: st_replacematrix("`tmp2'", st_matrix("`tmp2'"):/st_matrix("`B'")) } else if "`lnratio'"!="" { mata: st_replacematrix("`tmp2'", ln(st_matrix("`tmp2'"):/st_matrix("`B'"))) } else { matrix `tmp2' = `tmp2' - `B' } mat coleq `tmp2' = "`eq'" matrix `b' = nullmat(`b') , `tmp2' } matrix drop `B' `tmp' `tmp2' if "`se'"=="" { local eq0 if "`baseval'"!="+" local eq0 "`baseval'" foreach eq of local eqs { if "`baseval'"=="+" & "`eq0'"=="" local eq0 "`eq'" if "`eq'"=="`eq0'" { foreach eq1 of local eqs { if "`eq1'"=="`eq0'" continue mat `tmp' = nullmat(`tmp'), I(`bins') * -1 } } else { foreach eq1 of local eqs { if "`eq1'"=="`eq0'" continue if "`eq1'"=="`eq'" { mat `tmp' = nullmat(`tmp'), I(`bins') } else { mat `tmp' = nullmat(`tmp'), J(`bins', `bins', 0) } } } matrix `B' = nullmat(`B') \ `tmp' matrix drop `tmp' if "`baseval'"=="+" local eq0 "`eq'" } matrix `V' = `B'' * `V' * `B' matrix drop `B' } } // Post results if "`se'"=="" { local coln: colfullnames `b' mat coln `V' = `coln' mat rown `V' = `coln' } Erepost `b' `V' eret scalar k_eq = `k_eq' local title `"`e(title)'"' if "`ratio'"!="" local title: subinstr local title "Percentile" "Ratios in percentile" else if "`lnratio'"!="" local title: subinstr local title "Percentile" "Log. ratios in percentile" else local title: subinstr local title "Percentile" "Differences in percentile" eret local title "`title'" eret local contrast contrast eret local baseval "`baseval'" eret local ratio "`ratio'" eret local lnratio "`lnratio'" eret local _estimates_name "" end program Ratio_V, eclass tempname ecurrent _estimates hold `ecurrent', restore args ratio b V baseval if "`ratio'"=="lnratio" local ln ln local coln: colfullnames `b' mat coln `V' = `coln' mat rown `V' = `coln' local coefs: colnames `b' local coefs: list uniq coefs local eqs: coleq `b' local eqs: list uniq eqs eret post `b' `V' local coleq local coln local eq0 if "`baseval'"!="+" local eq0 "`baseval'" foreach eq of local eqs { if "`baseval'"=="+" & "`eq0'"=="" local eq0 "`eq'" if "`eq'"=="`eq0'" continue foreach c of local coefs { local nlcom `nlcom' (`ln'(_b[`eq':`c']/_b[`eq0':`c'])) local coleq `coleq' `eq' } local coln `coln' `coefs' if "`baseval'"=="+" local eq0 "`eq'" } qui nlcom `nlcom' mat `b' = r(b) mat coleq `b' = `coleq' mat coln `b' = `coln' mat `V' = r(V) end program Erepost, eclass args b V //backup existing e()'s tempvar sample gen byte `sample' = e(sample) local emacros: e(macros) foreach emacro of local emacros { local e_`emacro' `"`e(`emacro')'"' } local escalars: e(scalars) foreach escalar of local escalars { tempname e_`escalar' scalar `e_`escalar'' = e(`escalar') } local ematrices: e(matrices) local bV "b V" local ematrices: list ematrices - bV foreach ematrix of local ematrices { tempname e_`ematrix' matrix `e_`ematrix'' = e(`ematrix') } // post results eret post `b' `V', esample(`sample') foreach emacro of local emacros { eret local `emacro' `"`e_`emacro''"' } foreach escalar of local escalars { eret scalar `escalar' = scalar(`e_`escalar'') } foreach ematrix of local ematrices { eret matrix `ematrix' = `e_`ematrix'' } end version 11 mata mata set matastrict on mata: void pshare_stack_sort() // returns local eqs { real colvector S, p string colvector eqs if (st_local("sortgini")!="") { S = st_matrix(st_local("G")) } else { S = st_matrix("e(share)")' S = select(S, mod(1::rows(S),st_numscalar("e(bins)")):==0) } p = order(S, 1) if (st_local("descending")!="") { p = revorder(p) } eqs = tokens(st_local("eqs0"))' eqs = eqs[p] if (st_strscalar("e(over)")!="") { if (st_local("tlast")!="") { eqs = select(eqs, eqs:!="total") \ select(eqs, eqs:=="total") } else if (st_local("tfirst")!="") { eqs = select(eqs, eqs:=="total") \ select(eqs, eqs:!="total") } } st_local("eqs0", invtokens(eqs')) } void pshare_compute() // sets matrix S and local N, replaces uvars { string scalar rtype real scalar step, fw, se, N, ytot, wtot, wtotR, G, offset real colvector y, Y, w, W, p, P real matrix S pointer(real colvector) scalar pvar // setup rtype = st_local("rtype") step = (st_local("step")!="") fw = (st_local("weight")=="fweight") se = (st_local("uvars")!="") y = st_data(., st_local("y"), st_local("touse")) if (st_local("pvar")!="") { pvar = &(st_data(., st_local("pvar"), st_local("touse"))) } else pvar = &y if (st_local("weight")!="") { w = st_data(., st_local("exp"), st_local("touse")) Y = quadrunningsum(y:*w) W = quadrunningsum(w) N = rows(y) if (fw | (st_local("weight")=="iweight")) N = W[N] // sets N to sum of weights } else { w = 1 Y = quadrunningsum(y) N = rows(y) W = 1::N } if (st_local("normnum")!="") ytot = strtoreal(st_local("normnum")) else if (st_local("ytotal")!="") ytot = st_numscalar(st_local("ytotal")) else if (st_local("normvar")!="") ytot = quadsum(st_data(., st_local("normvar"), st_local("touse"))) else ytot = Y[rows(Y)] wtot = W[rows(W)] wtotR = (st_local("wtotal")!="" ? st_numscalar(st_local("wtotal")) : wtot) p = strtoreal(tokens(st_local("percentiles")))' P = (p/100) * wtot // compute shares S = (0 \ p) , (p \ 100) // ll, ul S = S, (S[,1] + (S[,2]-S[,1])/2), (S[,2]-S[,1])/100 // mid, prop S = S, J(rows(S), 3, .) // L_ll, L_ul, share if (se) S = S, J(rows(S), 2, .) if (step) { if (se) pshare_compute_step_se(S, *pvar, Y, W, P) else if (fw) pshare_compute_step_fw(S, Y, W, P) else S[|1,6 \ rows(S)-1,6|] = pshare_compute_step(Y, W, P)' } else { if (se) pshare_compute_ipolate_se(S, *pvar, Y, W, P) else S[|1,6 \ rows(S)-1,6|] = pshare_compute_ipolate(Y, W, P)' } S[rows(S), 6] = Y[rows(Y)] // scaling if (anyof(("density","proportion","percent"), rtype)) { S[,6] = S[,6] :/ ytot } if (rtype=="generalized") { S[,6] = S[,6] :/ wtot } S[,5] = 0 \ S[|1,6 \ rows(S)-1,6|] S[,7] = S[,6] - S[,5] if (rtype=="density") S[,7] = S[,7] :/ (S[,4] * (wtot/wtotR)) else if (rtype=="percent") S[,7] = S[,7] * 100 else if (rtype=="average") S[,7] = S[,7] :/ (wtot * S[,4]) // Gini if (st_local("gini")!="") { G = pshare_cindex((st_local("y0")!="" ? st_data(., st_local("y0"), st_local("touse")) : y), w, *pvar) } // set uvars if (se) { offset = strtoreal(st_local("offset")) S[,9] = S[,9] :+ offset pshare_compute_uvars(S, ytot, wtot, wtotR, rtype) S = S[|1,1 \ .,7|] if (st_local("over")!="") st_local("offset", strofreal(offset + rows(W), "%18.0g")) } // returns st_matrix(st_local("S"), S) st_local("N", strofreal(N, "%18.0g")) st_local("G", strofreal(G, "%18.0g")) } real rowvector pshare_compute_step( // determine lorenz ordinates from step function real colvector Y, real colvector W, real colvector P) { real scalar i, j, r real rowvector res r = rows(P) res = J(1, r, .) j = 1 for (i=1; i<=rows(Y); i++) { while (P[j]<=W[i]) { res[j] = Y[i] j++ if (j>r) break } if (j>r) break } return(res) } void pshare_compute_step_se( // determine lorenz ordinates from step function real matrix S, real colvector y, real colvector Y, real colvector W, real colvector P) { real scalar i, j, r r = rows(P) j = 1 for (i=1; i<=rows(Y); i++) { while (P[j]<=W[i]) { S[j, 6] = Y[i] S[j, 8] = y[i] // for SEs S[j, 9] = i // for SEs j++ if (j>r) break } if (j>r) break } } void pshare_compute_step_fw( // (interpolation to next integer) real matrix S, real colvector Y, real colvector W, real colvector P) { real scalar i, j, r, Y0, W0 r = rows(P) W0 = 0 Y0 = 0 j = 1 for (i=1; i<=rows(Y); i++) { while (P[j]<=W[i]) { if (i>1) { Y0 = Y[i-1] W0 = W[i-1] } S[j, 6] = Y0 + (Y[i] - Y0) * (ceil(P[j] - W0) / (W[i] - W0)) j++ if (j>r) break } if (j>r) break } } real rowvector pshare_compute_ipolate( // determine lorenz ordinates using interpolation real colvector Y, real colvector W, real colvector P) { real scalar i, j, r, Y0, W0 real rowvector res r = rows(P) res = J(1, r, .) W0 = Y0 = 0 // if (Y[1]<0) Y0 = Y[1] // no extrapolation if y can be negative j = 1 for (i=1; i<=rows(Y); i++) { while (P[j]<=W[i]) { if (i>1) { Y0 = Y[i-1] W0 = W[i-1] } res[j] = Y0 + (Y[i] - Y0) * ((P[j] - W0) / (W[i] - W0)) // note: W[i]>W0 always true, even of some weights are zero j++ if (j>r) break } if (j>r) break } return(res) } void pshare_compute_ipolate_se( // determine lorenz ordinates using interpolation real matrix S, real colvector y, real colvector Y, real colvector W, real colvector P) { real scalar i, j, r, y0, Y0, W0 r = rows(P) W0 = y0 = Y0 = 0 // if (Y[1]<0) y0 = Y0 = Y[1] // no extrapolation if y can be negative j = 1 for (i=1; i<=rows(Y); i++) { while (P[j]<=W[i]) { if (i>1) { y0 = y[i-1] Y0 = Y[i-1] W0 = W[i-1] } S[j, 6] = Y0 + (Y[i] - Y0) * ((P[j] - W0) / (W[i] - W0)) // note: W[i]>W0 always true, even if some weights are zero S[j, 8] = y0 + (y[i] - y0) * ((P[j] - W0) / (W[i] - W0)) // for SEs S[j, 9] = i // for SEs j++ if (j>r) break } if (j>r) break } } void pshare_compute_uvars( // sets the pseudo-variables for SE computation real matrix S, real scalar ytot, real scalar wtot, real scalar wtotR, string scalar rtype) { real scalar j real rowvector id, I0, I1 real colvector y, q, p, a, s, b, J real matrix uvars // estimate E[Y|Z = Q_p] in case of pvar() if (st_local("pvar")!="" & st_local("pvar")!=st_local("y0")) { pshare_compute_EY(S) // replaces S[,9] } // setup st_view(uvars, ., tokens(st_local("uvars")), st_local("touse0")) y = st_data(., (st_local("y0")!="" ? st_local("y0") : st_local("y")), st_local("touse0")) J = st_data(., st_local("touse"), st_local("touse0")) s = S[,7] if (anyof(("proportion", "percent", "density"), rtype)) { if (anyof(("proportion", "percent"), rtype)) { if (st_local("normnum")!="") { // 2 5 a = st_data(., st_local("avar"), st_local("touse0")) } else { // 0 1 3 4 a = (st_local("normvar")=="" ? y : // 0 3 st_data(., st_local("normvar"), st_local("touse0"))) // 1 4 if (st_local("normpop")=="") a = a :* J // 0 1 else { // 3 4 a = a :* st_data(., st_local("touseR"), st_local("touse0")) } } b = J(rows(S), 1, ytot) if (rtype=="percent") { s = s / 100 // better than changing a b = b / 100 } } else /*if (rtype=="density")*/ { if (st_local("normnum")!="") { // 2 5 a = st_data(., st_local("avar"), st_local("touse0")) if (st_local("normpop")!="") { // 5 a = (a :+ ((ytot/wtot) * J :- (ytot/wtotR) * st_data(., st_local("touseR"), st_local("touse0")))) * (wtot/wtotR) } } else { // 0 1 3 4 a = (st_local("normvar")=="" ? y : // 0 3 st_data(., st_local("normvar"), st_local("touse0"))) // 1 4 if (st_local("normpop")=="") a = a :* J // 0 1 else { // 3 4 a = (a :* st_data(., st_local("touseR"), st_local("touse0")) :+ ((ytot/wtot) * J :- (ytot/wtotR) * st_data(., st_local("touseR"), st_local("touse0")))) * (wtot/wtotR) } } s = s :* S[,4] b = (ytot * (wtot/wtotR)) :* S[,4] } } else if (rtype=="average") { a = J s = s :* S[,4] b = wtot * S[,4] } else if (rtype=="sum") { a = st_data(., st_local("avar"), st_local("touse0")) b = J(rows(S), 1, 1) } else /*if (rtype=="generalized")*/ { a = J b = J(rows(S), 1, wtot) } // fill in uvars id = 1::rows(y) j = 1 q = S[,8] p = S[,2] / 100 I1 = (id :<= S[j,9]) uvars[,j] = (((y :- q[j]) :* I1 :+ (p[j] * q[j])) :* J :- (a :* s[j])) / b[j] for (j=2; j<=rows(S)-1; j++) { swap(I0, I1) I1 = (id :<= S[j,9]) uvars[,j] = (((y :- q[j]) :* I1 :- (y :- q[j-1]) :* I0 :+ (p[j] * q[j] - p[j-1] * q[j-1])) :* J :- (a :* s[j])) / b[j] } j = rows(S) swap(I0, I1) uvars[,j] = ((y :- (y :- q[j-1]) :* I0 :- (p[j-1] * q[j-1])) :* J :- (a :* s[j])) / b[j] } void pshare_compute_EY(real matrix S) // estimate E[Y|Z = Q_p] using local linear regression { real colvector q real scalar preserve, rc string scalar y, z, w, s, at, yhat y = st_local("y0") z = st_local("pvar") w = st_local("exp") s = st_local("touse") q = S[|1,8 \ rows(S)-1,8|] preserve = (st_nobs()=.) continue if (mi==1) y[|mi \ i-1|] = J(i-mi, 1, y[i]) // missings at start else y[|mi \ i-1|] = y[mi-1] :+ (y[i] - y[mi-1]) :* (x[|mi \ i-1|] :- x[mi-1]) :/ (x[i] - x[mi-1]) mi = . continue } if (y[i]>=.) mi = i } if (mi==1) return(y) // all missing if (mi<.) y[|mi \ r|] = J(i-mi, 1, y[mi-1]) // missings at end return(y) } void pshare_compute_EY_mean(real matrix S) { // insert mean of y if there are still missings real scalar m if (!hasmissing(S[|1,8 \ rows(S)-1,8|])) return if (st_local("exp")=="") m = mean(st_data(., st_local("y0"), st_local("touse"))) else m = mean(st_data(., st_local("y0"), st_local("touse")), st_data(., st_local("exp"), st_local("touse"))) S[|1,8 \ rows(S)-1,8|] = editmissing(S[|1,8 \ rows(S)-1,8|], m) } // concentration (Gini) index (data must be sorted by y) real scalar pshare_cindex(real colvector x, real colvector w, real colvector y) { real matrix mv if (rows(x)<1) return(.) mv = pshare_meanvar0((x, pshare_ranks(y, w)), w) return(mv[3,1] * 2 / mv[1,1]) } real matrix pshare_meanvar0(real matrix X, real colvector w) { real rowvector CP real rowvector means real scalar n CP = quadcross(w,0, X,1) n = cols(CP) means = CP[|1\n-1|] :/ CP[n] if (missing(means)) return(means \ J(cols(X),cols(X),.)) return(means \ crossdev(X,0,means, w, X,0,means) :/ CP[n]) } real colvector pshare_ranks(real colvector x, real colvector w) { // assumes sorted input real scalar i, r, lr real colvector ranks if (rows(x)==0) return(J(0,1,.)) if (rows(w)!=1) ranks = runningsum(w) else ranks = (1::rows(x)) * w // normalize ranks ranks = ranks / ranks[rows(x)] // lowest rank in case of ties for(i=rows(x)-1; i>=1; i--) { if (x[i]==x[i+1]) ranks[i] = ranks[i+1] } // compute midpoints r = ranks[1] ranks[1] = r/2 lr = r for(i=2; i<=rows(ranks); i++) { if (x[i]==x[i-1]) ranks[i] = ranks[i-1] else { r = ranks[i] ranks[i] = lr + (r-lr)/2 lr = r } } return(ranks) } end exit