*! NJC 1.0.1 27 May 2007 *! NJC 1.0.0 22 May 2007 program firstdigit, rclass byable(recall) version 9.1 syntax varlist(numeric) [if] [in] /// [, BY(varname) ALLobs MISSing PERcent ] quietly { // what observations to use if "`allobs'" != "" marksample touse, novarlist else marksample touse if "`by'" != "" & "`missing'" == "" markout `touse' `by', strok count if `touse' if r(N) == 0 error 2000 // variable(s) or group(s) local nv : word count `varlist' // by option if "`by'" != "" { if `nv' > 1 { di as err /// "by() cannot be combined with `nv' variables" exit 198 } tempname stub separate `varlist' if `touse', /// by(`by') g(`stub') veryshortlabel `missing' local varlist "`r(varlist)'" local nv : word count `varlist' } // display preparation and initiation local nam = 0 foreach v of local varlist { if "`by'" != "" { local nam = /// max(`nam', length(`"`: variable label `v''"')) } else local nam = max(`nam', length("`v'")) } local col = min(28, `nam') + 1 local txt /// " n chi-sq. P-value digit observed expected" noi di _n as txt "{col `col'}`txt'" local line = `col' - 1 + length("`txt'") noi di "{hline `line'}" tempvar thisuse gen byte `thisuse' = . local ofmt = cond("`percent'" != "", "%11.2f", "%11.0f") // loop over variables (-separate- results if -by()-) foreach v of local varlist { local which = /// cond("`by'" != "", `"`: variable label `v''"', "`v'") local which = abbrev(`"`which'"', `col' - 1) replace `thisuse' = `touse' & `v' < . mata : fd_work("`v'", "`thisuse'", "`percent'") noi di as text `"`which'{col `col'}"' "`sp'" /// as res %8.0f r(N) /// as res %10.2f r(chisq) /// as res %9.4f r(p) /// as res " 1" `ofmt' r(obs1) %11.2f r(exp1) return scalar obs1 = r(obs1) return scalar exp1 = r(exp1) noi forval i = 2/9 { di as res "{col `col'}`sp'{space 34}`i'" /// `ofmt' r(obs`i') %11.2f r(exp`i') return scalar obs`i' = r(obs`i') return scalar exp`i' = r(exp`i') } return scalar p = r(p) return scalar chisq = r(chisq) return scalar N = r(N) local --nv if `nv' noi di " " } } end mata : void fd_work(string scalar varname, string scalar tousename, string scalar percent) { real colvector y, obs, exp real scalar n, i, chisq string scalar name y = st_data(., varname, tousename) n = rows(y) exp = obs = J(9, 1, 0) y = strtoreal(substr(strofreal(y), 1, 1)) for (i = 1; i <= 9; i++) { obs[i] = colsum(y :== i) exp[i] = n * log10(1 + 1/i) name = "r(obs" + strofreal(i) + ")" st_numscalar(name, percent == "" ? obs[i] : 100 * obs[i] / n) name = "r(exp" + strofreal(i) + ")" st_numscalar(name, percent == "" ? exp[i] : 100 * log10(1 + 1/i)) } chisq = colsum(((obs - exp):^2) :/ exp) st_numscalar("r(p)", chi2tail(8, chisq)) st_numscalar("r(chisq)", chisq) st_numscalar("r(N)", n) } end