Repository: numpy/numpy
Pull request: BUG: Fixed maximum relative error reporting in assert_allclose
Fixed maximum relative error reporting in `assert_allclose`.

In cases where the two arrays have zeros at the same positions it will
no longer report nan as the max relative error.

Fixes #13801.

<!-- Please be sure you are following the instructions in the dev guidelines
http://www.numpy.org/devdocs/dev/development_workflow.html
-->

<!-- We'd appreciate it if your commit message is properly formatted
http://www.numpy.org/devdocs/dev/development_workflow.html#writing-the-commit-message
-->


Problem statement:
assert_allclose reports NaN as max relative difference if arrays have zeros at the same positions
<!-- Please describe the issue in detail here, and fill in the fields below -->
If two different arrays share a zero at the same position the max relative difference will be reported as NaN by ```assert_allclose```. This is because in ```assert_array_compare``` the error is simply calculated as ```max_rel_error = (error / abs(y)).max()``` which can result in 0/0 = NaN. 
This hides the actual biggest relative difference and makes it more difficult to solve the problem in your code that is triggering the assertion. 
### Reproducing code example:

```python
import numpy as np
x = np.array([0, 1])
y = np.array([0, 2])
np.testing.assert_allclose(x, y, rtol=0.001)
```

### Output:
```
AssertionError: 
Not equal to tolerance rtol=0.001, atol=0

Mismatch: 50%
Max absolute difference: 1
Max relative difference: nan
 x: array([0, 1])
 y: array([0, 2])
```

### Expected output:
```
AssertionError: 
Not equal to tolerance rtol=0.001, atol=0

Mismatch: 50%
Max absolute difference: 1
Max relative difference: 0.5
 x: array([0, 1])
 y: array([0, 2])
```

### Numpy/Python version information:
1.16.4

Hints from issue discussion:

Commit message: Fixed maximum relative error reporting in assert_allclose:
In cases where the two arrays have zeros at the same positions it will
no longer report nan as the max relative error

Patch to review:
diff --git a/numpy/testing/_private/utils.py b/numpy/testing/_private/utils.py
index ead5d264d000..1f0c6012d1c2 100644
--- a/numpy/testing/_private/utils.py
+++ b/numpy/testing/_private/utils.py
@@ -704,6 +704,7 @@ def assert_array_compare(comparison, x, y, err_msg='', verbose=True,
                          equal_inf=True):
     __tracebackhide__ = True  # Hide traceback for py.test
     from numpy.core import array, array2string, isnan, inf, bool_, errstate
+    from numpy import nanmax
 
     x = array(x, copy=False, subok=True)
     y = array(y, copy=False, subok=True)
@@ -821,7 +822,12 @@ def func_assert_same_pos(x, y, func=isnan, hasval='nan'):
 
                     # note: this definition of relative error matches that one
                     # used by assert_allclose (found in np.isclose)
-                    max_rel_error = (error / abs(y)).max()
+                    # NaNs are disregarded to ensure that we report correct max
+                    # relative error between two arrays with zeros in the same position
+                    with warnings.catch_warnings():
+                        # Ignore warning in case of both array are full of zeros
+                        warnings.filterwarnings('ignore', r'All-NaN (slice|axis) encountered')
+                        max_rel_error = nanmax(error / abs(y))
                     if error.dtype == 'object':
                         remarks.append('Max relative difference: '
                                         + str(max_rel_error))
diff --git a/numpy/testing/tests/test_utils.py b/numpy/testing/tests/test_utils.py
index bf60772d3e11..1acc1158a455 100644
--- a/numpy/testing/tests/test_utils.py
+++ b/numpy/testing/tests/test_utils.py
@@ -880,6 +880,16 @@ def test_equal_nan_default(self):
         assert_array_less(a, b)
         assert_allclose(a, b)
 
+    def test_report_max_relative_error(self):
+        a = np.array([0, 1])
+        b = np.array([0, 2])
+
+        with pytest.raises(AssertionError) as exc_info:
+            assert_allclose(a, b)
+        msg = str(exc_info.value)
+        assert_('Mismatch: 50%\nMax absolute difference: 1\n'
+                'Max relative difference: 0.5' in msg)
+
 
 class TestArrayAlmostEqualNulp(object):
 


Code review comments:
diff_hunk: @@ -821,7 +822,12 @@ def func_assert_same_pos(x, y, func=isnan, hasval='nan'):
 
                     # note: this definition of relative error matches that one
                     # used by assert_allclose (found in np.isclose)
-                    max_rel_error = (error / abs(y)).max()
+                    # NaNs are disregarded to ensure that we report correct max
+                    # relative error between two arrays with zeros in the same position
+                    with warnings.catch_warnings():
path: numpy/testing/_private/utils.py
text: @user1:
Fixing the error is nice. Do we have a way to avoid the warnings context here? If necessary by using numpy internals. Warning contexts are not thread-safe, and I think we had problems before with such things because some projects used threading during testing (and if just because they test their own threading support).

@author:
We could use a stripped-down version of nanmax, but it would need to import ```_replace_nan``` from nanfunctions.py:
```python
rel_error = error / abs(y)
if type(rel_error) is ndarray and rel_error.dtype != object_:
    max_rel_error = fmax.reduce(rel_error, axis=None)
else:
    rel_error, mask = _replace_nan(rel_error, -inf)
    # If rel_error is full of NaNs print NaN as max relative error
    if mask is not None and all(mask):
        max_rel_error = array(nan)
    else:
        max_rel_error = amax(rel_error)
```
There don't seem to be any tests for ```assert_array_compare``` with object arrays though so I'm not sure if this solution is correct for them.

@user1:
@author maybe you could add such a test (it should be pretty small after all?). Sorry if this is a bit annoying. We could make a small private helper (reusing `_replace_nan`) with a comment that the only goal is to not have a warning maybe.

@author:
@user1 No problem! Is something like this ok?
```python
        def nanmax_nowarn(a, axis=None):
        """
        Same as nanmax, but doesn't raise a warning if the array is all NaN
        """
        from numpy import amax, fmax, all, object_
        from numpy.lib.nanfunctions import _replace_nan

        if type(a) is ndarray and a.dtype != object_:
            # Fast, but not safe for subclasses of ndarray, or object arrays,
            # which do not implement isnan (gh-9009), or fmax correctly (gh-8975)
            res = fmax.reduce(a, axis=axis)
        else:
            # Slow, but safe for subclasses of ndarray
            a, mask = _replace_nan(a, -inf)
            res = amax(a, axis=axis)
            if mask is None:
                return res

            # Check for all-NaN axis
            if mask.all():
                return array(nan)
        return res
```
One last question, hopefully: How should this handle subclasses of ```ndarray``` that don't implement ```all```? ```if mask.all()```would  raise ```NotImplementedError``` in cases like that.

@user1:
Those subclasses had issues before as well, no? The only would be duck types with `__array_function__` overrides, but those are brand new, and I not even sure if assert arrays compare doesn't convert to a numpy array by this time in any case.

I think this is right, the only thing I wonder if we can avoid the code duplication a bit....

@author:
> I think this is right, the only thing I wonder if we can avoid the code duplication a bit....

Well the simplest way to avoid rewriting ```nanmax``` would be to add a new parameter to it which allows the caller to decide whether they want warnings or not. I'm not sure if you'd want to change the user facing API just for this fix though.

@user1:
Yeah, I am very much unwilling to do that. Might be willing to make a private helper. Although the indiscretion looks a bit weird for the fmax branch. OTOH, that is exactly why... if that branching logic needs fixing up, it will not happen...

@user2:
@author @user1 is there any conclusion on a way forward here?

@author:
@user1 By "private helper" do you mean a helper called by nanmax?
It's kind of messy because of the branching inside nanmax, this is what I've come up with:
```python
def _nanmax(a, axis=None, out=None, keepdims=np._NoValue):
    """
    Helper function for nanmax.

    Returns
    -------
    nanmax : ndarray
        An array with the same shape as `a`, with the specified axis removed.
        If `a` is a 0-d array, or if axis is None, an ndarray scalar is
        returned.  The same dtype as `a` is returned.
    warning_msg : {string, None}
        None if no warning should be raised, the warning message if an all-NaN slice is encountered.
    """
    kwargs = {}
    if keepdims is not np._NoValue:
        kwargs['keepdims'] = keepdims
    warning_msg = None
    if type(a) is np.ndarray and a.dtype != np.object_:
        # Fast, but not safe for subclasses of ndarray, or object arrays,
        # which do not implement isnan (gh-9009), or fmax correctly (gh-8975)
        res = np.fmax.reduce(a, axis=axis, out=out, **kwargs)
        if np.isnan(res).any():
            warning_msg = "All-NaN slice encountered"
    else:
        # Slow, but safe for subclasses of ndarray
        a, mask = _replace_nan(a, -np.inf)
        res = np.amax(a, axis=axis, out=out, **kwargs)
        if mask is None:
            return res, None

        # Check for all-NaN axis
        mask = np.all(mask, axis=axis, **kwargs)
        if np.any(mask):
            res = _copyto(res, np.nan, mask)
            warning_msg = "All-NaN axis encountered"
    return res, warning_msg


@user3(_nanmax_dispatcher)
def nanmax(a, axis=None, out=None, keepdims=np._NoValue):
    """..."""
    res, warning_msg = _nanmax(a, axis, out, keepdims)
    if warning_msg is not None:
        warnings.warn(warning_msg, RuntimeWarning,
                      stacklevel=3)
    return res
```
Then inside ```assert_array_compare``` we can do this:
```python
rel_error = error / abs(y)
max_rel_error, _ = _nanmax(rel_error)
```

@user1:
Yeah, I think this is probably the best solution here, thanks for the followup.

@author:
Great! Should I change ```nanmin``` similarly for the sake of consistency or is that unnecessary?

@user1:
Hmmm, in principle there is a slight change in that we do not dispatch the max part anymore. I am willing to assume that is fine, also since most calls inside are ufuncs.

I do not think we need to do it for nanmin. Maybe add to the documentation that this helper only exists for use in `testing`.