# Can the empirical Hessian of an M-estimator be indefinite?

Jeffrey Wooldridge in his Econometric Analysis of Cross Section and Panel Data (page 357) says that the empirical Hessian "is not guaranteed to be positive definite, or even positive semidefinite, for the particular sample we are working with.".

This seems wrong to me as (numerical problems apart) the Hessian must be positive semidefinite as a result of the definition of the M-estimator as the value of the parameter which minimizes the objective function for the given sample and the well-known fact that at a (local) minimum the Hessian is positive semidefinite.

Is my argument right?

[EDIT: The statement has been removed in the 2nd ed. of the book. See comment.]

BACKGROUND Suppose that $\widehat \theta_N$ is an estimator obtained by minimizing $${1 \over N}\sum_{i=1}^N q(w_i,\theta),$$ where $w_i$ denotes the $i$-th observation.

Let's denote the Hessian of $q$ by $H$, $$H(q,\theta)_{ij}=\frac{\partial^2 q}{\partial \theta_i \partial \theta_j}$$

The asymptotic covariance of $\widehat \theta_n$ involves $E[H(q,\theta_0)]$ where $\theta_0$ is the true parameter value. One way to estimate it is to use the empirical Hesssian

$$\widehat H=\frac{1}{N}\sum_{i=1}^N H(w_i,\widehat \theta_n)$$

It is the definiteness of $\widehat H$ which is in question.

• @Jyotirmoy, what if the minimum happens at the boundary of your parameter space? – cardinal Feb 16 '11 at 20:18
• @cardinal. Your are right, my argument won't work in that case. But Wooldridge is considering the case where the minimum is in the interior. Isn't he wrong in that case? – Jyotirmoy Bhattacharya Feb 17 '11 at 4:02
• @Jyotirmoy, it can certainly be only positive semidefinite. Think of linear functions or a function where the set of minimum points forms a convex polytope. For a simpler example, consider any polynomial $f(x)=x^{2n}$ at $x = 0$. – cardinal Feb 17 '11 at 4:13
• @cardinal. True. What is troubling me is the phrase "even positive semidefinite" in the quoted statement. – Jyotirmoy Bhattacharya Feb 17 '11 at 5:04
• @Jyotirmoy, is there a specific form of the M-estimator given in the book that you could provide? Also give the parameter space under consideration. Maybe then we can figure out what the author had in mind. In general, I think we've already established that the author's assertion is correct. Placing further constraints on the form of $q$ or the parameter space being considered might alter that. – cardinal Feb 17 '11 at 14:47

I think you're right. Let's distill your argument to its essence:

1. $\widehat \theta_N$ minimizes the function $Q$ defined as $Q(\theta) = {1 \over N}\sum_{i=1}^N q(w_i,\theta).$

2. Let $H$ be the Hessian of $Q$, whence $H(\theta) = \frac{\partial^2 Q}{\partial \theta_i \partial \theta_j}$ by definition and this in turn, by linearity of differentiation, equals $\frac{1}{N}\sum_{i=1}^N H(w_i, \theta_n)$.

3. Assuming $\widehat \theta_N$ lies in the interior of the domain of $Q$, then $H(\widehat \theta_N)$ must be positive semi-definite.

This is merely a statement about the function $Q$: how it is defined is merely a distraction, except insofar as the assumed second order differentiability of $q$ with respect to its second argument ($\theta$) assures the second order differentiability of $Q$.

Finding M-estimators can be tricky. Consider these data provided by @mpiktas:

{1.168042, 0.3998378}, {1.807516, 0.5939584}, {1.384942, 3.6700205}, {1.327734, -3.3390724}, {1.602101, 4.1317608}, {1.604394, -1.9045958}, {1.124633, -3.0865249}, {1.294601, -1.8331763},{1.577610, 1.0865977}, { 1.630979, 0.7869717}


The R procedure to find the M-estimator with $q((x,y),\theta)=(y-c_1x^{c_2})^4$ produced the solution $(c_1, c_2)$ = $(-114.91316, -32.54386)$. The value of the objective function (the average of the $q$'s) at this point equals 62.3542. Here is a plot of the fit:

Here is a plot of the (log) objective function in a neighborhood of this fit:

Something is fishy here: the parameters of the fit are extremely far from the parameters used to simulate the data (near $(0.3, 0.2)$) and we do not seem to be at a minimum: we are in an extremely shallow valley that is sloping towards larger values of both parameters:

The negative determinant of the Hessian at this point confirms that this is not a local minimum! Nevertheless, when you look at the z-axis labels, you can see that this function is flat to five-digit precision within the entire region, because it equals a constant 4.1329 (the logarithm of 62.354). This probably led the R function minimizer (with its default tolerances) to conclude it was near a minimum.

In fact, the solution is far from this point. To be sure of finding it, I employed the computationally expensive but highly effective "Principal Axis" method in Mathematica, using 50-digit precision (base 10) to avoid possible numerical problems. It finds a minimum near $(c_1, c_2) = (0.02506, 7.55973)$ where the objective function has the value 58.292655: about 6% smaller than the "minimum" found by R. This minimum occurs in an extremely flat-looking section, but I can make it look (just barely) like a true minimum, with elliptical contours, by exaggerating the $c_2$ direction in the plot:

The contours range from 58.29266 in the middle all the way up to 58.29284 in the corners(!). Here's the 3D view (again of the log objective):

Here the Hessian is positive-definite: its eigenvalues are 55062.02 and 0.430978. Thus this point is a local minimum (and likely a global minimum). Here is the fit it corresponds to:

I think it's better than the other one. The parameter values are certainly more realistic and it's clear we're not going to be able to do much better with this family of curves.

There are useful lessons we can draw from this example:

1. Numerical optimization can be difficult, especially with nonlinear fitting and non-quadratic loss functions. Therefore:
2. Double-check results in as many ways as possible, including:
3. Graph the objective function whenever you can.
4. When numerical results appear to violate mathematical theorems, be extremely suspicious.
5. When statistical results are surprising--such as the surprising parameter values returned by the R code--be extra suspicious.
• +1, nice analysis. I think that is why Wooldridge included the remark. I still think it is possible to think of some example where the hessian will be indefinite. Artificially restricting the parameter space for example. In this example the parameter space is whole plane, that is why the local minimum will give semi-positive hessian. I think the time has come to write a nice email to Wooldridge to get his take on the question:) – mpiktas Feb 25 '11 at 18:00
• @mpiktas Yes, I'm sure there exist problems where an interior global minimum has an indefinite Hessian, yet where all parameters are identifiable. But it simply is not possible for the Hessian at a sufficiently smooth interior global minimum to be indefinite. This sort of thing has been proven again and again, such as in Milnor's Topology from a Differentiable Viewpoint. I suspect Wooldridge may have been misled by errant numerical "solutions." (The typos on the quoted page suggest it was written hastily, by the way.) – whuber Feb 25 '11 at 18:31
• even at the boundary, hessian will be positive? I'll check out the book, I see that I really lack extensive knowledge in this area. Classical theorems are very simple, so I assumed that there should not be something else very complicated. That maybe one of the reasons why I had so much difficulty answering the question. – mpiktas Feb 25 '11 at 19:19
• @mpiktas At the boundary the Hessian won't necessarily even be defined. The idea is this: if the Jacobian/Hessian/second derivative matrix is defined at a critical point, then in a neighborhood the function acts like the quadratic form determined by this matrix. If the matrix has positive and negative eigenvalues, the function must increase in some directions and decrease in others: it cannot be a local extremum. This is what concerned @Jyotirmoy about the quotation, which seems to contradict this basic property. – whuber Feb 25 '11 at 20:51
• Thanks you both you and @mpiktas for the very nice analysis. I would tend to agree with you that Wooldridge is confusing numerical difficulties with theoretical properties of the estimator. Let's see if there are any other answers. – Jyotirmoy Bhattacharya Feb 26 '11 at 4:22

The quotation in full can be found here. The estimate $\hat{\theta}_N$ is the solution of minimization problem (page 344):

\begin{align} \min_{\theta\in \Theta}N^{-1}\sum_{i=1}^Nq(w_i,\theta) \end{align}

If the solution $\hat{\theta}_N$ is interior point of $\Theta$, objective function is twice differentiable and gradient of the objective function is zero, then Hessian of the objective function (which is $\hat{H}$) is positive semi-definite.

Now what Wooldridge is saying that for given sample the empirical Hessian is not guaranteed to be positive definite or even positive semidefinite. This is true, since Wooldridge does not require that objective function $N^{-1}\sum_{i=1}^Nq(w_i,\theta)$ has nice properties, he requires that there exists a unique solution $\theta_0$ for

$$\min_{\theta\in\Theta}Eq(w,\theta).$$

So for given sample objective function $N^{-1}\sum_{i=1}^Nq(w_i,\theta)$ may be minimized on the boundary point of $\Theta$ in which Hessian of objective function needs not to be positive definite.

Further in his book Wooldridge gives an examples of estimates of Hessian which are guaranteed to be numerically positive definite. In practice non-positive definiteness of Hessian should indicate that solution is either on the boundary point or the algorithm failed to find the solution. Which usually is a further indication that the model fitted may be inappropriate for a given data.

Here is the numerical example. I generate non-linear least squares problem:

$$y_i=c_1x_i^{c_2}+\varepsilon_i$$

I take $X$ uniformly distributed in interval $[1,2]$ and $\varepsilon$ normal with zero mean and variance $\sigma^2$. I generated a sample of size 10, in R 2.11.1 using set.seed(3). Here is the link to the values of $x_i$ and $y_i$.

I chose the objective function square of usual non-linear least squares objective function:

$$q(w,\theta)=(y-c_1x_i^{c_2})^4$$

Here is the code in R for optimising function, its gradient and hessian.

##First set-up the epxressions for optimising function, its gradient and hessian.
##I use symbolic derivation of R to guard against human error
mt <- expression((y-c1*x^c2)^4)

##Evaluate the expressions on data to get the empirical values.
##Note there was a bug in previous version of the answer res should not be squared.
optf <- function(p) {
res <- eval(mt,list(y=y,x=x,c1=p[1],c2=p[2]))
mean(res)
}

gf <- function(p) {
evl <- list(y=y,x=x,c1=p[1],c2=p[2])
apply(res,2,mean)
}

hesf <- function(p) {
evl <- list(y=y,x=x,c1=p[1],c2=p[2])
res1 <- lapply(hessmt,function(l)sapply(l,function(ll)eval(ll,evl)))
res <- sapply(res1,function(l)apply(l,2,mean))
res
}


set.seed(3)
x <- runif(10,1,2)
y <- 0.3*x^0.2

> optf(c(0.3,0.2))
[1] 0
> gf(c(0.3,0.2))
[1] 0 0
> hesf(c(0.3,0.2))
[,1] [,2]
[1,]    0    0
[2,]    0    0
> eigen(hesf(c(0.3,0.2)))$values [1] 0 0  The hessian is zero, so it is positive semi-definite. Now for the values of$x$and$y$given in the link we get > df <- read.csv("badhessian.csv") > df x y 1 1.168042 0.3998378 2 1.807516 0.5939584 3 1.384942 3.6700205 4 1.327734 -3.3390724 5 1.602101 4.1317608 6 1.604394 -1.9045958 7 1.124633 -3.0865249 8 1.294601 -1.8331763 9 1.577610 1.0865977 10 1.630979 0.7869717 > x <- df$x
> y <- df$y > opt <- optim(c(1,1),optf,gr=gf,method="BFGS") > opt$par
[1] -114.91316  -32.54386
> gf(opt$par) [1] -0.0005795979 -0.0002399711 > hesf(opt$par)
[,1]         [,2]
[1,]  0.0002514806 -0.003670634
[2,] -0.0036706345  0.050998404
> eigen(hesf(opt$par))$values
[1]  5.126253e-02 -1.264959e-05


Gradient is zero, but the hessian is non positive.

Note: This is my third attempt to give an answer. I hope I finally managed to give precise mathematical statements, which eluded me in the previous versions.

• @mpiktas, That's some interesting notation there (I know it's not yours). A $w$ on the left-hand side and $y$ and $x$ on the right-hand side. I'm guessing $w = (x,y)$ or something like that. Also, I'm assuming the squaring should be happening to $y - m(x,\theta)$ and not just to $m(x,\theta)$. No? – cardinal Feb 18 '11 at 13:17
• @mpiktas, I'm not quite sure how to interpret your first sentence due to the wording. I can see two ways, one that I'd call correct and the other I wouldn't. Also, strictly speaking, I don't agree with the second sentence in your first paragraph. As I've shown above, it is possible to be at a local minimum in the interior of the parameter space without the Hessian being positive definite. – cardinal Feb 18 '11 at 13:24
• @cardinal, yes you are right. Wooldridge uses $w$ for consistency reasons, $y$ and $x$ is reserved for response and predictors throughout the book. In this example $w=(x,y)$. – mpiktas Feb 18 '11 at 13:26
• @cardinal, I fixed my wording. Now it should be ok. Thanks for pointing out the problem. – mpiktas Feb 18 '11 at 13:42
• @mptikas. Neither Wooldridge nor I are claiming that the Hessian has to be positive definite everywhere. My claim is that for an interior maximum the empirical Hessian has to be positive semidefinite as a necessary condition of a smooth function reaching its maximum. Wooldridge seems to be saying something different. – Jyotirmoy Bhattacharya Feb 24 '11 at 10:07

The hessian is indefinite at a saddle point. It’s possible that this may be the only stationary point in the interior of the parameter space.

Update: Let me elaborate. First, let’s assume that the empirical Hessian exists everywhere.

If $\hat{\theta}_n$ is a local (or even global) minimum of $\sum_i q(w_i, \cdot)$ and in the interior of the parameter space (assumed to be an open set) then necessarily the Hessian $(1/N) \sum_i H(w_i, \hat{\theta}_n)$ is positive semidefinite. If not, then $\hat{\theta}_n$ is not a local minimum. This follows from second order optimality conditions — locally $\sum_i q(w_i, \cdot)$ must not decrease in any directions away from $\hat{\theta}_n$.

One source of the confusion might the "working" definition of an M-estimator. Although in principle an M-estimator should be defined as $\arg\min_\theta \sum_i q(w_i, \theta)$, it might also be defined as a solution to the equation $$0 = \sum_i \dot{q}(w_i, \theta)\,,$$ where $\dot{q}$ is the gradient of $q(w, \theta)$ with respect to $\theta$. This is sometimes called the $\Psi$-type. In the latter case a solution of that equation need not be a local minimum. It can be a saddle point and in this case the Hessian would be indefinite.

Practically speaking, even a positive definite Hessian that is nearly singular or ill-conditioned would suggest that the estimator is poor and you have more to worry about than estimating its variance.

• could you adapt your answer so that it matches the notation of the question? To what is $x^2-y^2$ referring? Where does this get inserted into the equations given in the question? – probabilityislogic Feb 25 '11 at 8:44
• +1 Good points in the update, especially the last paragraph. When the Hessian is available--as is implicitly assumed throughout this discussion--one would automatically use its positive-definiteness as one of the criteria for testing any critical point and therefore this issue simply could not arise. This leads me to believe the Wooldridge quotation must concern the Hessian at a putative global minimum, not at a mere critical point. – whuber Feb 26 '11 at 13:52

There's been a lot of beating around the bush in this thread regarding whether the Hessian has to be positive (semi)definite at a local minimum. So I will make a clear statement on that.

Presuming the objective function and all constraint functions are twice continuously differentiable, then at any local minimum, the Hessian of the Lagrangian projected into the null space of the Jacobian of active constraints must be positive semidefinite. I.e., if $Z$ is a basis for the null space of the Jacobian of active constraints, then $Z^T*(\text{Hessian of Lagrangian})*Z$ must be positive semidefinite. This must be positive definite for a strict local minimum.

So the Hessian of the objective function in a constrained problem having active constraint(s) need not be positive semidefinite if there are active constraints.

Notes:

1) Active constraints consist of all equality constraints, plus inequality constraints which are satisfied with equality.

2) See the definition of the Lagrangian at https://www.encyclopediaofmath.org/index.php/Karush-Kuhn-Tucker_conditions .

3) If all constraints are linear, then the Hessian of the Lagrangian = Hessian of the objective function because the 2nd derivatives of linear functions are zero. But you still need to do the projection jazz if any of these constraints are active. Note that lower or upper bound constraints are particular cases of linear inequality constraints. If the only constraints which are active are bound constraints, the projection of the Hessian into the null space of the Jacobian of active constraints amounts to eliminating the rows and columns of the Hessian corresponding to those components on their bounds.

4) Because Lagrange multipliers of inactive constraints are zero, if there are no active constraints, the Hessian of the Lagrangian = the Hessian of the objective function, and the Identity matrix is a basis for the null space of the Jacobian of active constraints, which results in the simplification of the criterion being the familiar condition that the Hessian of the objective function be positive semidefinite at a local minimum (positive definite if a strict local minimum).

The positive answers above are true but they leave out the crucial identification assumption - if your model is not identified (or if it is only set identified) you might indeed, as Wooldridge correctly indicated, find yourself with a non-PSD empirical Hessian. Just run some non-toy psychometric / econometric model and see for yourself.

• Because this does not seem mathematically possible, could you offer a simple, clear example to demonstrate how the Hessian of a continuously twice-differentiable objective function could possibly fail to be PSD at a global minimum? – whuber Feb 12 '19 at 14:17