estat instructions show statistics after estimation. Many of those statistics are diagnostics or checks used to guage mannequin specification. Some statistics can be found in any case estimation instructions; others are command particular.
I illustrate how estat instructions work after which present how you can write a command-specific estat command for the mypoisson command that I’ve been creating.
That is the twenty eighth publish within the collection Programming an estimation command in Stata. I like to recommend that you simply begin initially. See Programming an estimation command in Stata: A map to posted entries for a map to all of the posts on this collection.
estat by instance
I wish to know whether or not the poisson mannequin in instance 1 is an inexpensive specification of the true course of.
Instance 1: A Poisson mannequin
. use accident3
. poisson accidents cvalue tickets site visitors
Iteration 0: log chance = -119.40895
Iteration 1: log chance = -118.11766
Iteration 2: log chance = -118.11672
Iteration 3: log chance = -118.11672
Poisson regression Variety of obs = 505
LR chi2(3) = 1215.60
Prob > chi2 = 0.0000
Log chance = -118.11672 Pseudo R2 = 0.8373
------------------------------------------------------------------------------
accidents | Coef. Std. Err. z P>|z| [95% Conf. Interval]
-------------+----------------------------------------------------------------
cvalue | .1936988 .1054813 1.84 0.066 -.0130408 .4004385
tickets | 1.528282 .0741663 20.61 0.000 1.382919 1.673646
site visitors | .0822459 .0259017 3.18 0.001 .0314795 .1330123
_cons | -8.838649 .5160757 -17.13 0.000 -9.850139 -7.82716
------------------------------------------------------------------------------
After any estimation command, estat summarize summarizes the variables used within the beforehand estimated mannequin.
Instance 2: estat summarize
. estat summarize
Estimation pattern poisson Variety of obs = 505
-------------------------------------------------------------------
Variable | Imply Std. Dev. Min Max
-------------+-----------------------------------------------------
accidents | .4712871 2.172992 0 20
cvalue | 2.467327 1.101686 1 4
tickets | 1.708911 2.022792 0 7
site visitors | 6.776931 2.390108 .3951611 9.998234
-------------------------------------------------------------------
estat gof studies a goodness-of-fit check after some estimation instructions. Instance 3 reveals that estat gof doesn’t reject the specification of the beforehand estimated mannequin.
Instance 3: estat gof
. estat gof
Deviance goodness-of-fit = 70.08114
Prob > chi2(501) = 1.0000
Pearson goodness-of-fit = 66.3997
Prob > chi2(501) = 1.0000
linktest implements one other specification check. linktest estimates the parameters of a mannequin, together with the linear prediction and its sq.. Below the null speculation of right specification, the coefficient on the sq. of the linear prediction ought to be zero. Instance 4 reveals that linktest rejects the beforehand match poisson specification.
Instance 4: linktest
. linktest
Iteration 0: log chance = -1050.1837
Iteration 1: log chance = -1004.8894 (backed up)
Iteration 2: log chance = -530.12448
Iteration 3: log chance = -128.25646
Iteration 4: log chance = -115.0814
Iteration 5: log chance = -111.70136
Iteration 6: log chance = -111.57302
Iteration 7: log chance = -111.57255
Iteration 8: log chance = -111.57255
Poisson regression Variety of obs = 505
LR chi2(2) = 1228.69
Prob > chi2 = 0.0000
Log chance = -111.57255 Pseudo R2 = 0.8463
------------------------------------------------------------------------------
accidents | Coef. Std. Err. z P>|z| [95% Conf. Interval]
-------------+----------------------------------------------------------------
_hat | 1.269188 .1134938 11.18 0.000 1.046744 1.491632
_hatsq | -.1135629 .036681 -3.10 0.002 -.1854563 -.0416695
_cons | .0873643 .1141074 0.77 0.444 -.1362821 .3110108
------------------------------------------------------------------------------
Instance 5 clarifies what linktest does.
Instance 5: What linktest does
. predict double xb, xb
. generate double xb2 = xb^2
. poisson accidents xb xb2
Iteration 0: log chance = -1050.1837
Iteration 1: log chance = -1004.8894 (backed up)
Iteration 2: log chance = -530.12475
Iteration 3: log chance = -128.25648
Iteration 4: log chance = -115.08141
Iteration 5: log chance = -111.70136
Iteration 6: log chance = -111.57302
Iteration 7: log chance = -111.57255
Iteration 8: log chance = -111.57255
Poisson regression Variety of obs = 505
LR chi2(2) = 1228.69
Prob > chi2 = 0.0000
Log chance = -111.57255 Pseudo R2 = 0.8463
------------------------------------------------------------------------------
accidents | Coef. Std. Err. z P>|z| [95% Conf. Interval]
-------------+----------------------------------------------------------------
xb | 1.269188 .1134938 11.18 0.000 1.046744 1.491632
xb2 | -.1135629 .036681 -3.10 0.002 -.1854563 -.0416695
_cons | .0873643 .1141074 0.77 0.444 -.1362821 .3110108
------------------------------------------------------------------------------
Within the subsequent part, I present how estat works after estimation instructions, and I write the command-specific estat command estat mylinktest, which performs a hyperlink check after mypoisson6, an up to date model of mypoisson5.
Hyperlink check outcomes utilizing mypossion5
mypoisson5.ado works with predict by calling mypoisson5_p.ado, as I mentioned in Programming an estimation command in Stata: Making predict work. In instance 6, I take advantage of mypoisson5 and its predict command to compute the predictions and run the Poisson regression wanted for a hyperlink check.
Instance 6: Hyperlink check items from mypoisson5
. mypoisson5 accidents cvalue tickets site visitors
Iteration 0: f(p) = -838.34841
Iteration 1: f(p) = -419.976
Iteration 2: f(p) = -145.89693
Iteration 3: f(p) = -121.07379
Iteration 4: f(p) = -118.12037
Iteration 5: f(p) = -118.11672
Iteration 6: f(p) = -118.11672
------------------------------------------------------------------------------
| Coef. Std. Err. z P>|z| [95% Conf. Interval]
-------------+----------------------------------------------------------------
cvalue | .1936988 .1054813 1.84 0.066 -.0130408 .4004385
tickets | 1.528282 .0741663 20.61 0.000 1.382919 1.673646
site visitors | .0822459 .0259017 3.18 0.001 .0314795 .1330123
_cons | -8.838649 .5160757 -17.13 0.000 -9.850139 -7.82716
------------------------------------------------------------------------------
. predict double myxb , xb
. generate double myxb2 = myxb^2
. mypoisson5 accidents myxb myxb2
Iteration 0: f(p) = -1003.1259
Iteration 1: f(p) = -641.52691
Iteration 2: f(p) = -185.06148
Iteration 3: f(p) = -137.97727
Iteration 4: f(p) = -121.82821
Iteration 5: f(p) = -114.13233
Iteration 6: f(p) = -111.62662
Iteration 7: f(p) = -111.57257
Iteration 8: f(p) = -111.57255
Iteration 9: f(p) = -111.57255
------------------------------------------------------------------------------
| Coef. Std. Err. z P>|z| [95% Conf. Interval]
-------------+----------------------------------------------------------------
myxb | 1.269188 .1134938 11.18 0.000 1.046744 1.491632
myxb2 | -.1135629 .036681 -3.10 0.002 -.1854563 -.0416695
_cons | .0873643 .1141074 0.77 0.444 -.1362821 .3110108
------------------------------------------------------------------------------
mypossion6 and estat_mylinktest
estat works very similar to predict: you retailer the title of your command-specific estat command within the e() outcomes of your estimation command. Apart from the title change, solely traces 44 and 51 in mypoisson6.ado differ from their counterparts in mypoisson5.ado.
*! model 6.0.0 18Oct2016
program outline mypoisson6, eclass sortpreserve
model 14.2
syntax varlist(numeric ts fv min=2) [if] [in] [, noCONStant vce(string) ]
marksample touse
_vce_parse `touse' , optlist(Strong) argoptlist(CLuster) : , vce(`vce')
native vce "`r(vce)'"
native clustervar "`r(cluster)'"
if "`vce'" == "strong" | "`vce'" == "cluster" {
native vcetype "Strong"
}
if "`clustervar'" != "" {
seize verify 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)
}
kind `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 publish `b' `V', esample(`touse') buildfvinfo depname(`depvar')
ereturn scalar N = `N'
ereturn scalar rank = `rank'
ereturn native vce "`vce'"
ereturn native vcetype "`vcetype'"
ereturn native clustvar "`clustervar'"
ereturn native predict "mypoisson6_p"
ereturn native estat_cmd "mypoisson6_estat"
ereturn native cmd "mypoisson6"
ereturn show
finish
program getcinfo, rclass
syntax varlist(ts fv), [ noCONStant ]
_rmcoll `varlist' , `fixed' develop
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, "gf2")
optimize_init_params(S, J(1, p, .01))
optimize_init_constraints(S, Ct)
b = optimize(S)
if (vcetype == "strong") {
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, mu
xb = X*b'
mu = exp(xb)
val = (-mu + y:*xb - lnfactorial(y))
if (todo>=1) {
grad = (y - mu):*X
}
if (todo==2) {
hess = -quadcross(X, mu, X)
}
}
finish
Line 44 now specifies choice depvar() on ereturn publish to retailer the title of the dependent variable in e(depvar). Line 51 is new; it shops mypoisson6_estat in e(estat_cmd) in order that estat can name the command-specific estat command mypoisson6_estat listed in code block 2.
*! model 1.0.0 18Oct2016
program mypoisson6_estat, rclass
model 14.2
if "`e(cmd)'" != "mypoisson6" {
error 301
}
gettoken subcmd relaxation : 0, parse(" ,")
if "`subcmd'"=="mylinktest" {
tempname eresults
tempvar xb xb2
native depvar = e(depvar)
predict double `xb' , xb
generate double `xb2' = `xb'^2
nobreak {
_estimates maintain `eresults'
mylinktestwork `depvar' `xb' `xb2'
native chi2 = e(lt_chi2)
native df = e(lt_df)
native p = e(lt_p)
_estimates unhold `eresults'
}
return clear
return scalar chi2 = `chi2'
return scalar df = `df'
return scalar p = `p'
}
else {
estat_default `0'
return add
}
finish
program mylinktestwork, eclass
syntax varlist(min=3 max=3)
tempvar b V
quietly mypoisson6 `varlist'
matrix `b' = e(b)
matrix colnames `b' = _hat _hatsq _cons
ereturn repost b = `b' , rename
ereturn show
quietly check _hatsq
native chi2 = r(chi2)
native df = r(df)
native p = r(p)
ereturn scalar lt_chi2 = `chi2'
ereturn scalar lt_df = `df'
ereturn scalar lt_p = `p'
finish
Traces 5–7 make sure that estat can solely name mypoisson6_estat when the lively e() outcomes had been produced by mypoisson6.
estat passes every little thing after the title “estat” to the command-specific estat command. For instance, if I sort
estat mylinktest
estat will go “mylinktest” to mypoisson6_estat.
If I sort
estat summarize, noheader
estat will go “summarize, noheader” to mypoisson6_estat. Line 9 makes use of gettoken to parse what estat passes; it would put the title of the estat subcommand within the native macro subcmd and every little thing that follows within the native macro relaxation.
Line 10 specifies that traces 11–28 will likely be executed if the consumer sorts estat mylinktest and that traces 31–32 will in any other case be executed. Traces 11–28 implement the hyperlink check. Traces 31–32 use estat_default to provide the outcomes for the not-command-specific estat instructions, like estat summarize.
Traces 11–16 put the required predictions into the non permanent variables xb and xb2. Line 17 ensures that traces 18–23 will likely be executed even when the consumer presses break whereas these traces are being executed. Line 18 shops the mypoisson6 e() leads to reminiscence, and line 23 places them again in e(), overwriting the e() outcomes utilized by mylinktestwork to compute and report the hyperlink check. Placing traces 18–23 in a nobreak block ensures that the unique mypoisson6 e()
outcomes are left lively, though mylinktestwork produces its personal e() outcomes.
Traces 36–53 implement the work routine mylinktestwork. Line 43 places the purpose estimates computed by working mypoisson6 on the dependent variable, the linear prediction, and its sq. within the Stata matrix saved within the native macro b. (For simplicity, let me name this vector b.) At this level, the column names of b are the names of the non permanent variables contained within the native macros `xb’ and `xb2′. These non permanent names usually are not informative. Line 43 places informative column names names on b. Line 44 makes use of ereturn repost with choice rename to retailer this new b in e(b) and to rename the row and column names of e(V) to be the column names in b. Line 45 shows the properly named link-test outcomes. Traces 46–52 compute and retailer the link-test outcomes which are retrieved and subsequently saved in r() by traces 20–22 and contours 25–28, respectively.
Line 15 makes use of predict, which is now dealt with by mypoisson6_p.ado.
*! model 1.0.0 18Oct2016
program outline mypoisson6_p
model 14.2
syntax newvarname [if] [in] , [ xb n ]
marksample touse, novarlist
native nopts : phrase rely `xb' `n'
if `nopts' >1 {
show "{err}just one statistic could also be specified"
exit 498
}
if `nopts' == 0 {
native n n
show "anticipated counts"
}
if "`xb'" != "" {
_predict `typlist' `varlist' if `touse' , xb
}
else {
tempvar xbv
quietly _predict double `xbv' if `touse' , xb
generate `typlist' `varlist' = exp(`xbv') if `touse'
}
finish
Instance 7 illustrates the method and reveals the output produced.
Instance 7: estat mylinktest
. mypoisson6 accidents cvalue tickets site visitors
Iteration 0: f(p) = -838.34841
Iteration 1: f(p) = -419.976
Iteration 2: f(p) = -145.89693
Iteration 3: f(p) = -121.07379
Iteration 4: f(p) = -118.12037
Iteration 5: f(p) = -118.11672
Iteration 6: f(p) = -118.11672
------------------------------------------------------------------------------
accidents | Coef. Std. Err. z P>|z| [95% Conf. Interval]
-------------+----------------------------------------------------------------
cvalue | .1936988 .1054813 1.84 0.066 -.0130408 .4004385
tickets | 1.528282 .0741663 20.61 0.000 1.382919 1.673646
site visitors | .0822459 .0259017 3.18 0.001 .0314795 .1330123
_cons | -8.838649 .5160757 -17.13 0.000 -9.850139 -7.82716
------------------------------------------------------------------------------
. estat mylinktest
------------------------------------------------------------------------------
accidents | Coef. Std. Err. z P>|z| [95% Conf. Interval]
-------------+----------------------------------------------------------------
_hat | 1.269188 .1134938 11.18 0.000 1.046744 1.491632
_hatsq | -.1135629 .036681 -3.10 0.002 -.1854563 -.0416695
_cons | .0873643 .1141074 0.77 0.444 -.1362821 .3110108
------------------------------------------------------------------------------
Instance 8 illustrates that mypoisson6_estat permits estat summarize to work.
Instance 8: estat summarize
. estat summarize
Estimation pattern mypoisson6 Variety of obs = 505
-------------------------------------------------------------------
Variable | Imply Std. Dev. Min Max
-------------+-----------------------------------------------------
accidents | .4712871 2.172992 0 20
cvalue | 2.467327 1.101686 1 4
tickets | 1.708911 2.022792 0 7
site visitors | 6.776931 2.390108 .3951611 9.998234
-------------------------------------------------------------------
Finished and undone
I illustrated how estat instructions work and how you can write a command-specific estat command for the mypoisson command that I’ve been creating.
