mypoisson3.ado provides choices for a strong or a cluster–sturdy estimator of the variance–covariance of the estimator (VCE) to mypoisson2.ado, which I mentioned in Programming an estimation command in Stata: Dealing with issue variables in a poisson command utilizing Mata. mypoisson3.ado parses the vce() choice utilizing the strategies I mentioned in Programming an estimation command in Stata: Including sturdy and cluster–sturdy VCEs to our Mata primarily based OLS command. Under, I present easy methods to use optimize() to compute the sturdy or cluster–sturdy VCE.
I solely talk about what’s new within the code for mypoisson3.ado, assuming that you’re accustomed to mypoisson2.ado.
That is the twenty-second submit within the sequence Programming an estimation command in Stata. I like to recommend that you just begin at first. See Programming an estimation command in Stata: A map to posted entries for a map to all of the posts on this sequence.
A poisson command with choices for a strong or a cluster–sturdy VCE
mypoisson3 computes Poisson-regression leads to Mata. The syntax of the mypoisson3 command is
mypoisson3 depvar indepvars [if] [in] [, vce(robust | cluster clustervar) noconstant]
the place indepvars can include issue variables or time-series variables.
Within the the rest of this submit, I talk about the code for mypoisson3.ado. I like to recommend that you just click on on the filename to obtain the code. To keep away from scrolling, view the code within the Do-file Editor, or your favourite textual content editor, to see the road numbers.
*! model 3.0.0 21Feb2016
program outline mypoisson3, eclass sortpreserve
model 14
syntax varlist(numeric ts fv min=2) [if] [in] [, noCONStant vce(string) ]
marksample touse
_vce_parse `touse' , optlist(Sturdy) argoptlist(CLuster) : , vce(`vce')
native vce "`r(vce)'"
native clustervar "`r(cluster)'"
if "`vce'" == "sturdy" | "`vce'" == "cluster" {
native vcetype "Sturdy"
}
if "`clustervar'" != "" {
seize affirm numeric variable `clustervar'
if _rc {
show in crimson "invalid vce() choice"
show in crimson "cluster variable {bf:`clustervar'} is " ///
"string variable as a substitute of a numeric variable"
exit(198)
}
type `clustervar'
}
gettoken depvar indepvars : varlist
_fv_check_depvar `depvar'
tempname b mo V N rank
getcinfo `indepvars' , `fixed'
native cnames "`r(cnames)'"
matrix `mo' = r(mo)
mata: mywork("`depvar'", "`cnames'", "`touse'", "`fixed'", ///
"`b'", "`V'", "`N'", "`rank'", "`mo'", "`vce'", "`clustervar'")
if "`fixed'" == "" {
native cnames "`cnames' _cons"
}
matrix colnames `b' = `cnames'
matrix colnames `V' = `cnames'
matrix rownames `V' = `cnames'
ereturn submit `b' `V', esample(`touse') buildfvinfo
ereturn scalar N = `N'
ereturn scalar rank = `rank'
ereturn native vce "`vce'"
ereturn native vcetype "`vcetype'"
ereturn native clustvar "`clustervar'"
ereturn native cmd "mypoisson3"
ereturn show
finish
program getcinfo, rclass
syntax varlist(ts fv), [ noCONStant ]
_rmcoll `varlist' , `fixed' increase
native cnames `r(varlist)'
native p : phrase rely `cnames'
if "`fixed'" == "" {
native p = `p' + 1
native cons _cons
}
tempname b mo
matrix `b' = J(1, `p', 0)
matrix colnames `b' = `cnames' `cons'
_ms_omit_info `b'
matrix `mo' = r(omit)
return native cnames "`cnames'"
return matrix mo = `mo'
finish
mata:
void mywork( string scalar depvar, string scalar indepvars,
string scalar touse, string scalar fixed,
string scalar bname, string scalar Vname,
string scalar nname, string scalar rname,
string scalar mo,
string scalar vcetype, string scalar clustervar)
{
actual vector y, b
actual matrix X, V, Ct
actual scalar n, p, rank
y = st_data(., depvar, touse)
n = rows(y)
X = st_data(., indepvars, touse)
if (fixed == "") {
X = X,J(n, 1, 1)
}
p = cols(X)
Ct = makeCt(mo)
S = optimize_init()
optimize_init_argument(S, 1, y)
optimize_init_argument(S, 2, X)
optimize_init_evaluator(S, &plleval3())
optimize_init_evaluatortype(S, "gf0")
optimize_init_params(S, J(1, p, .01))
optimize_init_constraints(S, Ct)
b = optimize(S)
if (vcetype == "sturdy") {
V = optimize_result_V_robust(S)
}
else if (vcetype == "cluster") {
cvar = st_data(., clustervar, touse)
optimize_init_cluster(S, cvar)
V = optimize_result_V_robust(S)
}
else { // vcetype should IID
V = optimize_result_V_oim(S)
}
rank = p - diag0cnt(invsym(V))
st_matrix(bname, b)
st_matrix(Vname, V)
st_numscalar(nname, n)
st_numscalar(rname, rank)
}
actual matrix makeCt(string scalar mo)
{
actual vector mo_v
actual scalar ko, j, p
mo_v = st_matrix(mo)
p = cols(mo_v)
ko = sum(mo_v)
if (ko>0) {
Ct = J(0, p, .)
for(j=1; j<=p; j++) {
if (mo_v[j]==1) {
Ct = Ct e(j, p)
}
}
Ct = Ct, J(ko, 1, 0)
}
else {
Ct = J(0,p+1,.)
}
return(Ct)
}
void plleval3(actual scalar todo, actual vector b, ///
actual vector y, actual matrix X, ///
val, grad, hess)
{
actual vector xb
xb = X*b'
val = (-exp(xb) + y:*xb - lnfactorial(y))
}
finish
Only some traces of mypoisson3.ado differ from their counterparts in mypoisson2.ado, and I put these adjustments into 4 teams.
- Line 5 permits vce() on the syntax command, and features 8–23 parse this selection.
I mentioned the strategies utilized in these adjustments in Programming an estimation com-
mand in Stata: Including sturdy and cluster–sturdy VCEs to our Mata primarily based OLS command, once I used them in myregress12.ado. These traces- put the required VCE within the native macro vce;
- put a label for the required VCE within the native macro vcetype;
- put the title of a specified cluster variable within the native macro clustervar, and
- deal with any errors when the person misspecifies the vce() choice.
- Line 35 passes the contents of the native macros vce and clustervar to the Mata work operate mywork().
- Strains 47–49 retailer the native macros vce, vcetype, and clustvar in e() outcomes.
- Line 84 parses the brand new arguments vcetype and clustervar. The string scalar vcetype comprises the kind of VCE to be estimated, and the string scalar clustervar comprises the title of the Stata variable containing the clusters, if specified.
- Strains 112–122 use the contents of vcetype to return an OIM, a strong, or a cluster–sturdy estimator of the VCE.
The contents of vcetype decide which optimize() operate is known as to compute the estimated VCE. If vcetype comprises sturdy, line 113 makes use of optimize_result_V_robust() to compute a strong estimator of the VCE. If vcetype comprises cluster, traces 116 and 117 put a replica of the Stata cluster variable within the optimize object, after which line 118 makes use of optimize_result_V_robust() to compute a cluster–sturdy estimator of the VCE. Lastly, if vcetype is empty, line 121 makes use of optimize_result_V_oim() to compute the default correct-specification estimator of the VCE.
The output in examples 1 and a couple of confirms that mypoisson3 produces the identical outcomes as poisson when the choice vce(cluster id) is specified.
Instance 1: mypoisson3 outcomes
. clear all
. use accident3
. mypoisson3 accidents cvalue i.children visitors, vce(cluster id)
Iteration 0: f(p) = -847.19028
Iteration 1: f(p) = -573.7331
Iteration 2: f(p) = -545.76673
Iteration 3: f(p) = -545.11357
Iteration 4: f(p) = -545.10898
Iteration 5: f(p) = -545.10898
(Std. Err. adjusted for clustering on id)
------------------------------------------------------------------------------
| Sturdy
| Coef. Std. Err. z P>|z| [95% Conf. Interval]
-------------+----------------------------------------------------------------
cvalue | -.6582924 .1128794 -5.83 0.000 -.8795319 -.4370529
|
children |
1 | -1.662351 .4309205 -3.86 0.000 -2.50694 -.8177623
2 | -1.574691 .4164515 -3.78 0.000 -2.390921 -.7584611
3 | -3.233933 .4685643 -6.90 0.000 -4.152302 -2.315564
|
visitors | .1383976 .0876168 1.58 0.114 -.0333282 .3101235
_cons | .7157579 .5970943 1.20 0.231 -.4545254 1.886041
------------------------------------------------------------------------------
Instance 2: poisson outcomes
. poisson accidents cvalue i.children visitors, vce(cluster id)
Iteration 0: log pseudolikelihood = -546.35782
Iteration 1: log pseudolikelihood = -545.11016
Iteration 2: log pseudolikelihood = -545.10898
Iteration 3: log pseudolikelihood = -545.10898
Poisson regression Variety of obs = 505
Wald chi2(5) = 118.06
Prob > chi2 = 0.0000
Log pseudolikelihood = -545.10898 Pseudo R2 = 0.2491
(Std. Err. adjusted for 285 clusters in id)
------------------------------------------------------------------------------
| Sturdy
accidents | Coef. Std. Err. z P>|z| [95% Conf. Interval]
-------------+----------------------------------------------------------------
cvalue | -.6582924 .1128793 -5.83 0.000 -.8795317 -.437053
|
children |
1 | -1.662351 .4309205 -3.86 0.000 -2.50694 -.8177622
2 | -1.574691 .4164515 -3.78 0.000 -2.390921 -.758461
3 | -3.233932 .4685642 -6.90 0.000 -4.152301 -2.315563
|
visitors | .1383977 .0876167 1.58 0.114 -.0333279 .3101232
_cons | .7157576 .597093 1.20 0.231 -.4545232 1.886038
------------------------------------------------------------------------------
Performed and undone
I mentioned mypoisson3, which has choices for a strong or a cluster–sturdy estimator of the variance–covariance of the estimator. In my subsequent submit, I talk about easy methods to have the evaluator operate compute the derivatives to hurry up the optimization.
