 |
[Date Prev][Date Next][Thread Prev][Thread Next][Date index][Thread index]
Re: st: condition of Hessian
Teresa Dale Nelson <[email protected]> asked for an example of how to
compute the condition number of the Hessian from an estimator reporting an
observed information matrix (OIM) variance.
Here is an example, using a silly Poisson regression:
. set mem 200m
(204800k)
. webuse nlswork
(National Longitudinal Survey. Young Women 14-26 years of age in 1968)
. poisson grade nev_mar not_smsa union
Iteration 0: log likelihood = -46328.9
Iteration 1: log likelihood = -46328.9
Poisson regression Number of obs = 19222
LR chi2(3) = 168.52
Prob > chi2 = 0.0000
Log likelihood = -46328.9 Pseudo R2 = 0.0018
------------------------------------------------------------------------------
grade | Coef. Std. Err. z P>|z| [95% Conf. Interval]
-------------+----------------------------------------------------------------
nev_mar | .0162206 .0050936 3.18 0.001 .0062374 .0262038
not_smsa | -.0507212 .0045589 -11.13 0.000 -.0596564 -.0417859
union | .0209936 .0047465 4.42 0.000 .0116907 .0302966
_cons | 2.552624 .0028711 889.06 0.000 2.546997 2.558251
------------------------------------------------------------------------------
.
. mata: V1 = st_matrix("e(V)")
. mata: cond(V1)
8.835376744
This condition number is close to one, as far as condition numbers go,
indicating a well-defined problem.
This condition number is based on the 2 norm,
cond = = norm(A, 2) * norm(A^(-1), 2)
where norm(A, 2) is the largest singular value of A. This calculation will
yield the same answer for A and A^(-1), ignoring precision problems.
Here are some illustrative computations
. mata: norm(V1,2)
.0000300734
. mata: svdsv(V1)
1
+---------------+
1 | .0000300734 |
2 | .0000244259 |
3 | .0000195973 |
4 | 3.40375e-06 |
+---------------+
. mata: norm(V1,2)*norm(cholinv(V1),2)
8.835376744
Now I claimed that 8.835 was close to one. Here I construct a badly
conditioned problem and calculate the condition number of the OIM variance
estimator.
. generate prob = union + .001*uniform()
(9296 missing values generated)
. poisson grade nev_mar not_smsa union prob
Iteration 0: log likelihood = -46328.894
Iteration 1: log likelihood = -46328.894
Poisson regression Number of obs = 19222
LR chi2(4) = 168.53
Prob > chi2 = 0.0000
Log likelihood = -46328.894 Pseudo R2 = 0.0018
------------------------------------------------------------------------------
grade | Coef. Std. Err. z P>|z| [95% Conf. Interval]
-------------+----------------------------------------------------------------
nev_mar | .0162228 .0050936 3.18 0.001 .0062395 .0262061
not_smsa | -.0507209 .0045589 -11.13 0.000 -.0596561 -.0417857
union | -.7705006 6.98085 -0.11 0.912 -14.45272 12.91171
prob | .7914966 6.980869 0.11 0.910 -12.89076 14.47375
_cons | 2.552229 .0045136 565.46 0.000 2.543383 2.561076
------------------------------------------------------------------------------
. mata: V2 = st_matrix("e(V)")
. mata: cond(V2)
30757045.59
This condition number is not close to one, indicating that the problem is
not well defined.
The data must be on the same scale to use condition numbers as a metric to
measure identification.
Here I illustrate by rescaling the union indicator variable. I compute the
correlation matrix before and after the rescaling to illustrate that the
correlations are invariant to the rescaling.
. corr grade nev_mar not_smsa union
(obs=19222)
| grade nev_mar not_smsa union
-------------+------------------------------------
grade | 1.0000
nev_mar | 0.0450 1.0000
not_smsa | -0.1275 -0.0759 1.0000
union | 0.0573 0.0195 -0.0684 1.0000
. replace union = 100000*union
union was byte now long
(4510 real changes made)
. corr grade nev_mar not_smsa union
(obs=19222)
| grade nev_mar not_smsa union
-------------+------------------------------------
grade | 1.0000
nev_mar | 0.0450 1.0000
not_smsa | -0.1275 -0.0759 1.0000
union | 0.0573 0.0195 -0.0684 1.0000
Now I refit the model and calculate the condition number of the OIM
variance.
.
. poisson grade nev_mar not_smsa union
Iteration 0: log likelihood = -46328.9
Iteration 1: log likelihood = -46328.9
Poisson regression Number of obs = 19222
LR chi2(3) = 168.52
Prob > chi2 = 0.0000
Log likelihood = -46328.9 Pseudo R2 = 0.0018
------------------------------------------------------------------------------
grade | Coef. Std. Err. z P>|z| [95% Conf. Interval]
-------------+----------------------------------------------------------------
nev_mar | .0162206 .0050936 3.18 0.001 .0062374 .0262038
not_smsa | -.0507212 .0045589 -11.13 0.000 -.0596564 -.0417859
union | 2.10e-07 4.75e-08 4.42 0.000 1.17e-07 3.03e-07
_cons | 2.552624 .0028711 889.06 0.000 2.546997 2.558251
------------------------------------------------------------------------------
.
. mata: V2 = st_matrix("e(V)")
. mata: cond(V2)
1.70881e+10
The rescaling has made this well-defined problem appear ill defined.
Drukker and Wiggins (2004) discuss these issues using real data and a
different norm. Drukker and Wiggins (2004) place the problem in terms of
parameter identification.
I hope that this helps.
--David
References
---------
Drukker, D. M. and V. Wiggins. 2004. "Verifying the Solution from a
Nonlinear Solver: A Case Study: Comment". The American Economic Review
94(1):397-399
*
* For searches and help try:
* http://www.stata.com/support/faqs/res/findit.html
* http://www.stata.com/support/statalist/faq
* http://www.ats.ucla.edu/stat/stata/