**Recent addition:** Inspired by DimaPasechnik's and Matt F.'s comments about sum of squares decompositions, I tried the following very natural idea: Try to find $f$ of the form $f(x,y)=A(x,y)^2+B(x,y)^2$, where $A$ and $B$ are relatively prime quadratic polynomials whose zero sets intersect in $4$ points, which then of course are local minima, and attempt to get a fifth local minimum. By a shift we may look for it in $(0,0)$. And indeed, there are examples which can be checked instantly by hand: One possibility is
\begin{align}
A &= x^2 - y - 6\\
B &= 8x^2 - y^2 - 3y + 2.
\end{align}
The curves $A=0$ and $B=0$ intersect in the four points
\begin{equation}
(\pm1, -5)\text{ and } (\pm4, 10).
\end{equation}
As
\begin{equation}
f(x, y) = 40 + 20x^2 + 6y^2 + \text{terms of degree }\ge3,
\end{equation}
we see that $(0, 0)$ is a fifth local minimum.

**Original answer:** The quartic
\begin{equation}
f = 3x^4 - 5x^2y^2 + 5y^4 + 20x^2y - 32y^3 + x^2 + 2y^2,
\end{equation}
which is even in $x$, has $5$ isolated local minima in
\begin{equation}
(0, 0), (\pm2, -1), (\pm2, 5).
\end{equation}
In order to prove the claim, one checks that the partial derivatives $\frac{\partial f}{\partial x}$ and $\frac{\partial f}{\partial y}$ vanish and that the Hessian is positive definite in these points, where the latter is equivalent to the traces and the determinants being positive.

If one carries out the calculation by hand, then the symmetry $f(x,y)=f(-x,y)$ requires only to check $3$ points. Alternatively, the following Sage code verifies the example:

```
R.<x, y> = QQ[]
f = 3*x^4-5*x^2*y^2+5*y^4+20*x^2*y-32*y^3+x^2+2*y^2
fx, fy = f.derivative(x), f.derivative(y)
fxx = fx.derivative(x)
fxy = fx.derivative(y)
fyy = fy.derivative(y)
tra, det = fxx+fyy, fxx*fyy-fxy^2
cps = ideal([fx, fy]).variety()
minima = [p for p in cps if tra.subs(p) > 0 and det.subs(p) > 0]
print(len(minima) == 5)
```

*Remark 1:* As pointed out by Jap88 in the comments (see also his nice visualizations at https://math.stackexchange.com/questions/4620663), previous examples are difficult to visualize. This holds true for this one too. The reason is that some of the saddles points tend to be very close to the minima. The red and blue lines are the curves given by $\frac{\partial f}{\partial x}=0$ and $\frac{\partial f}{\partial y}=0$, and the green dots are the minima. As the first picture doesn't show well the situation around $(0,0)$, the second one zooms into this area:

*Remark 2:* Another rather short example is
\begin{equation}
f = x^4 - 3x^2y^2 + 3y^4 - 2x^2y + 3y^3 + x^2 - 2y^2.
\end{equation}
Here, the $5$ local minima are irrational, so the exact algebraic verification by hand is a little more involved. Or, if one relies on Sage again, then one simply replaces the first two lines of the code above with

```
R.<x, y> = AA[]
f = x^4 - 3*x^2*y^2 + 3*y^4 - 2*x^2*y + 3*y^3 + x^2 - 2*y^2
```

Here AA is the (exact) field of real algebraic numbers. The image (notation as above) is

*Remark 3:* Conjecture 5.4 in this paper by Durfee et.al predicts that a quartic has at most $4$ isolated local minima. So this is a counterexample.

*Remark 4*: Theorem 3.1.6. of the paper Critical points of real polynomials ... by Shustin seems to claim that $5$ minima are possible. However, the "proof" given in Section 3.13 merely says "Assume that $d \le 4$. Then all the index distributions except [...] can be easily realized as it was explained above". I don't know what the author means by "above", as all the previous cases handle degrees $d\ge5$.