3
$\begingroup$

In working on an estimation problem, I created three Gaussian-like functions. The ML estimate could be seen as the argmax of the product of the three functions. In trying to get the code to run faster, I discovered complementary functions to ArgMax: NArgMax, and FindArgMax. However, I can not get the functions NArgMax and FindArgMaxto agree with the ArgMax solution. For reference, I believe that the ArgMax solution is correct, and the other two are not. Without pasting all of my code, here is a MWE that demonstrates the effect:

test1[x_, y_] := 1/.01^2*Exp[-(x^2 + y^2)/0.01^2]
test2[x_, y_] := 1/.1^2*Exp[-((x - 5)^2 + (y - 3)^2)/0.1^2]
ArgMax[test1[x, y] test2[x, y], {x, y}] // AbsoluteTiming
FindArgMax[
  test1[x, y]*test2[x, y], {{x, 2.5}, {y, 1.5}}] // AbsoluteTiming
NArgMax[test1[x, y] test2[x, y], {x, y}] // AbsoluteTiming

The results are

{0.0388864, {0.049505, 0.029703}}
{0.102196, {2.5, 1.5}}
{0.0687933, {1.02414*10^30, 8.98768*10^29}}

However, it is worth noting that FindArgMax throws the warning: FindArgMax::fmgz: Encountered a gradient that is effectively zero. The result returned may not be a maximum; it may be a minimum or a saddle point. >>. NArgMax does not produce any warnings or errors.

I believe this warning may be the reason that neither FindArgMax nor NArgMax find the proper solution: FindArgMax simply considers the initial estimate to be the final solution (as it warns), but I do not know why NArgMax returns such large values for its solution. However, ArgMax does correctly identify the maximum.

Any help in resolving this would be much appreciated.

||improve this question
$\endgroup$
  • $\begingroup$ How is AbsoluteTiming relevant to the question? $\endgroup$ – Taiki Nov 3 '15 at 3:36
  • $\begingroup$ It isn't really relevant to the question, it was more me trying to figure out which would be faster, if they were in agreement. I didn't remove it because AbsoluteTiming is well documented and I figured it didn't really matter. I can remove it, if you think it really detracts from the question. $\endgroup$ – Michael Witt Nov 3 '15 at 3:38
  • $\begingroup$ Your test function test1[x, y] * test2[x, y] produces ridiculously small numbers for any x and y. I guess none of the results are correct. $\endgroup$ – Taiki Nov 3 '15 at 3:44
  • $\begingroup$ No, I'm fairly certain that the result produced by ArgMax is correct. The size of the numeric solution shouldn't matter. If you let the variance of test1 equal the variance of test2, the ML estimate should be the midpoint of the two centroids, which is what it produces in that case; and FindArgMax also would 'find' that result, although that is its initial estimate. $\endgroup$ – Michael Witt Nov 3 '15 at 3:52
  • 2
    $\begingroup$ try using exact values (100^2 instead of 1/.01^2) and increasing workingprecision $\endgroup$ – george2079 Nov 3 '15 at 4:01
2
$\begingroup$

As pointed out by george2079 in the comment, it is adviseable to use arbitrary-precision numbers in this case, so we can set WorkingPrecision however high.

test1[x_, y_] := 1/(1/100)^2*Exp[-(x^2 + y^2)/(1/100)^2]
test2[x_, y_] := 1/(1/10)^2*Exp[-((x - 5)^2 + (y - 3)^2)/(1/10)^2]
f[x_, y_] := test1[x_, y_] * test2[x_, y_]

Then, for example,

NArgMax[f[x, y], {x, y}, WorkingPrecision -> 10]

gives a result that agrees with that given by ArgMax.

As to FindArgMax, it's a little more tricky. Even

FindArgMax[f[x, y], {x, y}, WorkingPrecision -> 10000]

doesn't help. We can improve that by manually supplying the gradient and Hessian:

grad = D[f[x, y], {{x, y}}]
hess = D[f[x, y], {{x, y}}, {{x, y}}]
FindArgMax[
  f[x, y],
  {x, y},
  Gradient -> grad, 
  Method -> {"Newton", Hessian -> hess},
  WorkingPrecision -> 10000
]

FindArgMax now gives a result that agrees with that given by ArgMax.

||improve this answer
$\endgroup$
  • $\begingroup$ That 'fix' for FindArgMax still requires ridiculously high WorkingPrecision. Simply changing the estimate is a simpler and easier way of getting FindArgMax to be in agreement. Furthermore, there is the issue that I pointed out, that even setting WorkingPrecision -> 2 for NArgMax allows it to find the solution properly. Why would that help? A precision of 2 digits is extremely low, yet it allows NArgMax to converge. $\endgroup$ – Michael Witt Nov 3 '15 at 4:42
  • $\begingroup$ I don't know why NArgMax works with low WorkingPrecision and FindArgMax requires high WorkingPrecision. What did you mean by saying that 'that fix ... still requires ridiculously high WorkingPrecision'? Does the fix not resolve your problem as asked for in the question? $\endgroup$ – Taiki Nov 3 '15 at 4:51
  • $\begingroup$ Well I mean, it works now. But the default WorkingPrecision is MachinePrecision, which is roughly 15.95. So setting WorkingPrecision to 2 shouldn't really help, yet it does. So, yes it resolves the problem, but only raises more questions. $\endgroup$ – Michael Witt Nov 3 '15 at 4:57
  • $\begingroup$ Indeed it's interesting to investigate further why that's the case. Also strange is that NArgMax[f[x, y], {x, y}, WorkingPrecision -> $MachinePrecision] gives the result that agrees with ArgMax, while NArgMax[f[x, y], {x, y}, WorkingPrecision -> MachinePrecision] doesn't! $\endgroup$ – Taiki Nov 3 '15 at 5:03

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.