diff --git a/tests/test_axis.py b/tests/test_axis.py
index 089694c2e..344080b4a 100644
--- a/tests/test_axis.py
+++ b/tests/test_axis.py
@@ -3,7 +3,6 @@
 
 import numpy as np
 import pytest
-from numpy.testing import assert_allclose, assert_array_equal
 from pytest import approx
 
 import boost_histogram as bh
@@ -247,8 +246,8 @@ def test_getitem(self):
         a = bh.axis.Regular(2, 1.0, 2.0)
         ref = [1.0, 1.5, 2.0]
         for i in range(2):
-            assert_allclose(a.bin(i), ref[i : i + 2])
-            assert_allclose(a[i], ref[i : i + 2])
+            assert a.bin(i) == approx(ref[i : i + 2])
+            assert a[i] == approx(ref[i : i + 2])
 
         assert a[-1] == a[1]
         with pytest.raises(IndexError):
@@ -257,8 +256,8 @@ def test_getitem(self):
         assert a.bin(-1)[0] == -np.inf
         assert a.bin(2)[1] == np.inf
 
-        assert_allclose(a[bh.underflow], a.bin(-1))
-        assert_allclose(a[bh.overflow], a.bin(2))
+        assert a[bh.underflow] == approx(a.bin(-1))
+        assert a[bh.overflow] == approx(a.bin(2))
 
         with pytest.raises(IndexError):
             a.bin(-2)
@@ -268,7 +267,7 @@ def test_getitem(self):
     def test_iter(self):
         a = bh.axis.Regular(2, 1.0, 2.0)
         ref = [(1.0, 1.5), (1.5, 2.0)]
-        assert_allclose(a, ref)
+        assert a == approx(ref)
 
     def test_index(self):
         a = bh.axis.Regular(4, 1.0, 2.0)
@@ -356,9 +355,9 @@ def test_pow_transform(self):
 
     def test_edges_centers_widths(self):
         a = bh.axis.Regular(2, 0, 1)
-        assert_allclose(a.edges, [0, 0.5, 1])
-        assert_allclose(a.centers, [0.25, 0.75])
-        assert_allclose(a.widths, [0.5, 0.5])
+        assert a.edges == approx([0, 0.5, 1])
+        assert a.centers == approx([0.25, 0.75])
+        assert a.widths == approx([0.5, 0.5])
 
 
 class TestCircular(Axis):
@@ -412,8 +411,8 @@ def test_getitem(self):
         a = bh.axis.Regular(2, 1, 1 + np.pi * 2, circular=True)
         ref = [1.0, 1.0 + np.pi, 1.0 + 2.0 * np.pi]
         for i in range(2):
-            assert_allclose(a.bin(i), ref[i : i + 2])
-            assert_allclose(a[i], ref[i : i + 2])
+            assert a.bin(i) == approx(ref[i : i + 2])
+            assert a[i] == approx(ref[i : i + 2])
 
         assert a[-1] == a[1]
         with pytest.raises(IndexError):
@@ -422,7 +421,7 @@ def test_getitem(self):
         with pytest.raises(IndexError):
             a[bh.underflow]
 
-        assert_allclose(a[bh.overflow], a.bin(2))
+        assert a[bh.overflow] == approx(a.bin(2))
 
         assert a.bin(2)[0] == approx(1 + 2 * np.pi)
         assert a.bin(2)[1] == approx(1 + 3 * np.pi)
@@ -435,7 +434,7 @@ def test_getitem(self):
     def test_iter(self):
         a = bh.axis.Regular(2, 1, 2, circular=True)
         ref = [(1, 1.5), (1.5, 2)]
-        assert_allclose(a, ref)
+        assert a == approx(ref)
 
     def test_index(self):
         a = bh.axis.Regular(4, 1, 1 + np.pi * 2, circular=True)
@@ -456,9 +455,9 @@ def test_index(self):
 
     def test_edges_centers_widths(self):
         a = bh.axis.Regular(2, 0, 1, circular=True)
-        assert_allclose(a.edges, [0, 0.5, 1])
-        assert_allclose(a.centers, [0.25, 0.75])
-        assert_allclose(a.widths, [0.5, 0.5])
+        assert a.edges == approx([0, 0.5, 1])
+        assert a.centers == approx([0.25, 0.75])
+        assert a.widths == approx([0.5, 0.5])
 
 
 class TestVariable(Axis):
@@ -526,15 +525,15 @@ def test_getitem(self):
         a = bh.axis.Variable(ref)
 
         for i in range(2):
-            assert_allclose(a.bin(i), ref[i : i + 2])
-            assert_allclose(a[i], ref[i : i + 2])
+            assert a.bin(i) == approx(ref[i : i + 2])
+            assert a[i] == approx(ref[i : i + 2])
 
         assert a[-1] == a[1]
         with pytest.raises(IndexError):
             a[2]
 
-        assert_allclose(a[bh.underflow], a.bin(-1))
-        assert_allclose(a[bh.overflow], a.bin(2))
+        assert a[bh.underflow] == approx(a.bin(-1))
+        assert a[bh.overflow] == approx(a.bin(2))
 
         assert a.bin(-1)[0] == -np.inf
         assert a.bin(-1)[1] == ref[0]
@@ -551,7 +550,7 @@ def test_iter(self):
         ref = [-0.1, 0.2, 0.3]
         a = bh.axis.Variable(ref)
         for i, bin in enumerate(a):
-            assert_array_equal(bin, ref[i : i + 2])
+            assert bin == approx(ref[i : i + 2])
 
     def test_index(self):
         a = bh.axis.Variable([-0.1, 0.2, 0.3])
@@ -569,9 +568,9 @@ def test_index(self):
 
     def test_edges_centers_widths(self):
         a = bh.axis.Variable([0, 1, 3])
-        assert_allclose(a.edges, [0, 1, 3])
-        assert_allclose(a.centers, [0.5, 2])
-        assert_allclose(a.widths, [1, 2])
+        assert a.edges == approx([0, 1, 3])
+        assert a.centers == approx([0.5, 2])
+        assert a.widths == approx([1, 2])
 
 
 class TestInteger:
@@ -667,13 +666,13 @@ def test_getitem(self):
         assert a.bin(-1) == -2
         assert a.bin(4) == 3
 
-        assert_allclose(a[bh.underflow], a.bin(-1))
-        assert_allclose(a[bh.overflow], a.bin(4))
+        assert a[bh.underflow] == approx(a.bin(-1))
+        assert a[bh.overflow] == approx(a.bin(4))
 
     def test_iter(self):
         a = bh.axis.Integer(-1, 3)
         ref = (-1, 0, 1, 2)
-        assert_array_equal(a, ref)
+        assert a == approx(ref)
 
     def test_index(self):
         a = bh.axis.Integer(-1, 3)
@@ -688,9 +687,9 @@ def test_index(self):
 
     def test_edges_centers_widths(self):
         a = bh.axis.Integer(1, 3)
-        assert_allclose(a.edges, [1, 2, 3])
-        assert_allclose(a.centers, [1.5, 2.5])
-        assert_allclose(a.widths, [1, 1])
+        assert a.edges == approx([1, 2, 3])
+        assert a.centers == approx([1.5, 2.5])
+        assert a.widths == approx([1, 1])
 
 
 class TestCategory(Axis):
@@ -793,7 +792,7 @@ def test_getitem(self, ref, growth):
     def test_iter(self, ref, growth):
         Cat = bh.axis.StrCategory if isinstance(ref[0], str) else bh.axis.IntCategory
         a = Cat(ref, growth=growth)
-        assert_array_equal(a, ref)
+        assert a == approx(ref)
 
     @pytest.mark.parametrize(
         "ref", ([1, 2, 3, 4], ("A", "B", "C", "D")), ids=("int", "str")
@@ -803,8 +802,8 @@ def test_index(self, ref, growth):
         a = Cat(ref, growth=growth)
         for i, r in enumerate(ref):
             assert a.index(r) == i
-        assert_array_equal(a.index(ref), [0, 1, 2, 3])
-        assert_array_equal(a.index(np.reshape(ref, (2, 2))), [[0, 1], [2, 3]])
+        assert a.index(ref) == approx([0, 1, 2, 3])
+        assert a.index(np.reshape(ref, (2, 2))) == approx([[0, 1], [2, 3]])
 
         if isinstance(ref[0], str):
             with pytest.raises(KeyError):
@@ -821,12 +820,10 @@ def test_value(self, ref, growth):
         a = Cat(ref, growth=growth)
         for i, r in enumerate(ref):
             assert a.value(i) == r
-        assert_array_equal(a.value(range(3)), ref)
+        assert a.value(range(3)) == approx(ref)
         assert a.value(3) is None
-        assert_array_equal(a.value((0, 3)), [ref[0], None])
-        assert_array_equal(
-            a.value(np.array((0, 1, 2, 3))), [ref[0], ref[1], ref[2], None]
-        )
+        assert a.value((0, 3)) == approx([ref[0], None])
+        assert a.value(np.array((0, 1, 2, 3))) == approx([ref[0], ref[1], ref[2], None])
         # may be added in the future
         with pytest.raises(ValueError):
             a.value([[2], [2]])
@@ -835,9 +832,9 @@ def test_value(self, ref, growth):
     def test_edges_centers_widths(self, ref, growth):
         Cat = bh.axis.StrCategory if isinstance(ref[0], str) else bh.axis.IntCategory
         a = Cat(ref, growth=growth)
-        assert_allclose(a.edges, [0, 1, 2, 3])
-        assert_allclose(a.centers, [0.5, 1.5, 2.5])
-        assert_allclose(a.widths, [1, 1, 1])
+        assert a.edges == approx([0, 1, 2, 3])
+        assert a.centers == approx([0.5, 1.5, 2.5])
+        assert a.widths == approx([1, 1, 1])
 
 
 class TestBoolean:
@@ -891,7 +888,7 @@ def test_getitem(self):
     def test_iter(self):
         a = bh.axis.Boolean()
         ref = (False, True)
-        assert_array_equal(a, ref)
+        assert a == approx(ref)
 
     def test_index(self):
         a = bh.axis.Boolean()
@@ -900,6 +897,6 @@ def test_index(self):
 
     def test_edges_centers_widths(self):
         a = bh.axis.Boolean()
-        assert_allclose(a.edges, [0.0, 1.0, 2.0])
-        assert_allclose(a.centers, [0.5, 1.5])
-        assert_allclose(a.widths, [1, 1])
+        assert a.edges == approx([0.0, 1.0, 2.0])
+        assert a.centers == approx([0.5, 1.5])
+        assert a.widths == approx([1, 1])
diff --git a/tests/test_benchmark_1d.py b/tests/test_benchmark_1d.py
index 9a59535b5..50eedb303 100644
--- a/tests/test_benchmark_1d.py
+++ b/tests/test_benchmark_1d.py
@@ -2,7 +2,7 @@
 
 import numpy as np
 import pytest
-from numpy.testing import assert_allclose, assert_array_equal
+from pytest import approx
 
 import boost_histogram as bh
 
@@ -31,7 +31,7 @@
 @pytest.mark.parametrize("dtype", vals)
 def test_numpy_1d(benchmark, dtype):
     result, _ = benchmark(np.histogram, vals[dtype], bins=bins, range=ranges)
-    assert_array_equal(result, answer[dtype])
+    assert result == approx(answer[dtype])
 
 
 def make_and_run_hist(flow, storage, vals):
@@ -47,4 +47,4 @@ def make_and_run_hist(flow, storage, vals):
 @pytest.mark.parametrize("storage", STORAGES)
 def test_boost_1d(benchmark, flow, storage, dtype):
     result = benchmark(make_and_run_hist, flow, storage, vals[dtype])
-    assert_allclose(result[:-1], answer[dtype][:-1], atol=2)
+    assert result[:-1] == approx(answer[dtype][:-1], atol=2)
diff --git a/tests/test_benchmark_2d.py b/tests/test_benchmark_2d.py
index 8db2ac197..1e58354b5 100644
--- a/tests/test_benchmark_2d.py
+++ b/tests/test_benchmark_2d.py
@@ -2,7 +2,7 @@
 
 import numpy as np
 import pytest
-from numpy.testing import assert_array_equal
+from pytest import approx
 
 import boost_histogram as bh
 
@@ -34,7 +34,7 @@
 @pytest.mark.parametrize("dtype", vals)
 def test_numpy_perf_2d(benchmark, dtype):
     result, _, _ = benchmark(np.histogram2d, *vals[dtype], bins=bins, range=ranges)
-    assert_array_equal(result, answer[dtype])
+    assert result == approx(answer[dtype])
 
 
 def make_and_run_hist(flow, storage, vals):
@@ -53,4 +53,4 @@ def make_and_run_hist(flow, storage, vals):
 @pytest.mark.parametrize("storage", STORAGES)
 def test_2d(benchmark, flow, storage, dtype):
     result = benchmark(make_and_run_hist, flow, storage, vals[dtype])
-    assert_array_equal(result[:-1, :-1], answer[dtype][:-1, :-1])
+    assert result[:-1, :-1] == approx(answer[dtype][:-1, :-1])
diff --git a/tests/test_histogram.py b/tests/test_histogram.py
index 678f0449a..b9c7adf3d 100644
--- a/tests/test_histogram.py
+++ b/tests/test_histogram.py
@@ -8,7 +8,6 @@
 
 import numpy as np
 import pytest
-from numpy.testing import assert_array_equal
 from pytest import approx
 
 import boost_histogram as bh
@@ -164,7 +163,9 @@ def test_setting(count_single_storage):
     assert h[9] == 5
     assert h[bh.overflow] == 6
 
-    assert_array_equal(h.view(flow=True), [1, 2, 3, 0, 0, 0, 4, 0, 0, 0, 5, 6])
+    assert np.asarray(h.view(flow=True)) == approx(
+        np.asarray([1, 2, 3, 0, 0, 0, 4, 0, 0, 0, 5, 6])
+    )
 
 
 def test_growth():
@@ -600,11 +601,11 @@ def test_shrink_1d():
     h = bh.Histogram(bh.axis.Regular(20, 1, 5))
     h.fill(1.1)
     hs = h[{0: slice(bh.loc(1), bh.loc(2))}]
-    assert_array_equal(hs.view(), [1, 0, 0, 0, 0])
+    assert np.asarray(hs.view()) == approx(np.asarray([1, 0, 0, 0, 0]))
 
     d = OrderedDict({0: slice(bh.loc(1), bh.loc(2))})
     hs = h[d]
-    assert_array_equal(hs.view(), [1, 0, 0, 0, 0])
+    assert np.asarray(hs.view()) == approx(np.asarray([1, 0, 0, 0, 0]))
 
 
 def test_rebin_1d():
@@ -612,17 +613,17 @@ def test_rebin_1d():
     h.fill(1.1)
 
     hs = h[{0: slice(None, None, bh.rebin(4))}]
-    assert_array_equal(hs.view(), [1, 0, 0, 0, 0])
+    assert np.asarray(hs.view()) == approx(np.asarray([1, 0, 0, 0, 0]))
 
     hs = h[{0: bh.rebin(4)}]
-    assert_array_equal(hs.view(), [1, 0, 0, 0, 0])
+    assert np.asarray(hs.view()) == approx(np.asarray([1, 0, 0, 0, 0]))
 
 
 def test_shrink_rebin_1d():
     h = bh.Histogram(bh.axis.Regular(20, 0, 4))
     h.fill(1.1)
     hs = h[{0: slice(bh.loc(1), bh.loc(3), bh.rebin(2))}]
-    assert_array_equal(hs.view(), [1, 0, 0, 0, 0])
+    assert np.asarray(hs.view()) == approx(np.asarray([1, 0, 0, 0, 0]))
 
 
 def test_rebin_nd():
@@ -718,7 +719,7 @@ def test_fill_bool_not_bool():
 
     h.fill([0, 1, 1, 7, -3])
 
-    assert_array_equal(h.view(), [1, 4])
+    assert np.asarray(h.view()) == approx(np.asarray([1, 4]))
 
 
 def test_pick_bool():
@@ -727,22 +728,22 @@ def test_pick_bool():
     h.fill([True, True, False, False], [True, False, True, True])
     h.fill([True, True, True], True)
 
-    assert_array_equal(h[True, :].view(), [1, 4])
-    assert_array_equal(h[False, :].view(), [0, 2])
-    assert_array_equal(h[:, False].view(), [0, 1])
-    assert_array_equal(h[:, True].view(), [2, 4])
+    assert np.asarray(h[True, :].view()) == approx(np.asarray([1, 4]))
+    assert np.asarray(h[False, :].view()) == approx(np.asarray([0, 2]))
+    assert np.asarray(h[:, False].view()) == approx(np.asarray([0, 1]))
+    assert np.asarray(h[:, True].view()) == approx(np.asarray([2, 4]))
 
 
 def test_slice_bool():
     h = bh.Histogram(bh.axis.Boolean())
     h.fill([0, 0, 0, 1, 3, 4, -2])
 
-    assert_array_equal(h.view(), [3, 4])
-    assert_array_equal(h[1:].view(), [4])
-    assert_array_equal(h[:1].view(), [3])
+    assert np.asarray(h.view()) == approx(np.asarray([3, 4]))
+    assert np.asarray(h[1:].view()) == approx(np.asarray([4]))
+    assert np.asarray(h[:1].view()) == approx(np.asarray([3]))
 
-    assert_array_equal(h[:1].axes[0].centers, [0.5])
-    assert_array_equal(h[1:].axes[0].centers, [1.5])
+    assert np.asarray(h[:1].axes[0].centers) == approx(np.asarray([0.5]))
+    assert np.asarray(h[1:].axes[0].centers) == approx(np.asarray([1.5]))
 
 
 def test_pickle_bool():
@@ -770,7 +771,7 @@ def test_pickle_bool():
     assert repr(a) == repr(b)
     assert str(a) == str(b)
     assert a == b
-    assert_array_equal(a.view(), b.view())
+    assert np.asarray(a.view()) == approx(np.asarray(b.view()))
 
 
 # NumPy tests
@@ -786,20 +787,20 @@ def test_numpy_conversion_0():
 
     for t in (c, v):
         assert t.dtype == np.double  # CLASSIC: np.uint8
-        assert_array_equal(t, (1, 5, 0))
+        assert t == approx((1, 5, 0))
 
     for _ in range(10):
         a.fill(2)
     # copy does not change, but view does
-    assert_array_equal(c, (1, 5, 0))
-    assert_array_equal(v, (1, 5, 10))
+    assert c == approx((1, 5, 0))
+    assert v == approx((1, 5, 10))
 
     for _ in range(255):
         a.fill(1)
     c = np.array(a)
 
     assert c.dtype == np.double  # CLASSIC: np.uint16
-    assert_array_equal(c, (1, 260, 10))
+    assert c == approx((1, 260, 10))
     # view does not follow underlying switch in word size
     # assert not np.all(c, v)
 
@@ -812,8 +813,8 @@ def test_numpy_conversion_1():
     c = np.array(h)  # a copy
     v = np.asarray(h)  # a view
     assert c.dtype == np.double  # CLASSIC: np.float64
-    assert_array_equal(c, np.array((0, 30, 0)))
-    assert_array_equal(v, c)
+    assert c == approx(np.array((0, 30, 0)))
+    assert v == approx(c)
 
 
 def test_numpy_conversion_2():
@@ -836,13 +837,13 @@ def test_numpy_conversion_2():
             for k in range(a.axes[2].extent):
                 d[i, j, k] = a[i, j, k]
 
-    assert_array_equal(d, r)
+    assert d == approx(r)
 
     c = np.array(a)  # a copy
     v = np.asarray(a)  # a view
 
-    assert_array_equal(c, r)
-    assert_array_equal(v, r)
+    assert c == approx(r)
+    assert v == approx(r)
 
 
 def test_numpy_conversion_3():
@@ -861,7 +862,7 @@ def test_numpy_conversion_3():
                 r[i, j, k] = i + j + k
     c = a.view(flow=True)
 
-    assert_array_equal(c, r)
+    assert c == approx(r)
 
     assert a.sum() == approx(144)
     assert a.sum(flow=True) == approx(720)
@@ -927,7 +928,7 @@ def ia(*args):
     a.fill(np.array(1))  # 0-dim arrays work
     a.fill(ia(-1, 0, 1, 2))
     a.fill((2, 1, 0, -1))
-    assert_array_equal(a.view(True), [2, 2, 3, 2])
+    assert np.asarray(a.view(True)) == approx(np.asarray([2, 2, 3, 2]))
 
     with pytest.raises(ValueError):
         a.fill(np.empty((2, 2)))
@@ -944,13 +945,13 @@ def ia(*args):
     b = bh.Histogram(bh.axis.Regular(3, 0, 3))
     b.fill(fa(0, 0, 1, 2))
     b.fill(ia(1, 0, 2, 2))
-    assert_array_equal(b.view(True), [0, 3, 2, 3, 0])
+    assert np.asarray(b.view(True)) == approx(np.asarray([0, 3, 2, 3, 0]))
 
     c = bh.Histogram(
         bh.axis.Integer(0, 2, underflow=False, overflow=False), bh.axis.Regular(2, 0, 2)
     )
     c.fill(ia(-1, 0, 1), fa(-1.0, 1.5, 0.5))
-    assert_array_equal(c.view(True), [[0, 0, 1, 0], [0, 1, 0, 0]])
+    assert np.asarray(c.view(True)) == approx(np.asarray([[0, 0, 1, 0], [0, 1, 0, 0]]))
     # we don't support: assert a[[1, 1]].value, 0
 
     with pytest.raises(ValueError):
@@ -960,11 +961,11 @@ def ia(*args):
 
     # this broadcasts
     c.fill([1, 0], -1)
-    assert_array_equal(c.view(True), [[1, 0, 1, 0], [1, 1, 0, 0]])
+    assert np.asarray(c.view(True)) == approx(np.asarray([[1, 0, 1, 0], [1, 1, 0, 0]]))
     c.fill([1, 0], 0)
-    assert_array_equal(c.view(True), [[1, 1, 1, 0], [1, 2, 0, 0]])
+    assert np.asarray(c.view(True)) == approx(np.asarray([[1, 1, 1, 0], [1, 2, 0, 0]]))
     c.fill(0, [-1, 0.5, 1.5, 2.5])
-    assert_array_equal(c.view(True), [[2, 2, 2, 1], [1, 2, 0, 0]])
+    assert np.asarray(c.view(True)) == approx(np.asarray([[2, 2, 2, 1], [1, 2, 0, 0]]))
 
     with pytest.raises(IndexError):
         c[1]
@@ -1046,10 +1047,10 @@ def test_fill_with_sequence_2():
     a.fill("A")
     a.fill(np.array("B"))  # 0-dim array is also accepted
     a.fill(("A", "B", "C"))
-    assert_array_equal(a.view(True), [2, 2, 1])
+    assert np.asarray(a.view(True)) == approx(np.asarray([2, 2, 1]))
     a.fill(np.array(("D", "B", "A"), dtype="S5"))
     a.fill(np.array(("D", "B", "A"), dtype="U1"))
-    assert_array_equal(a.view(True), [4, 4, 3])
+    assert np.asarray(a.view(True)) == approx(np.asarray([4, 4, 3]))
 
     with pytest.raises(ValueError):
         a.fill(np.array((("B", "A"), ("C", "A"))))  # ndim == 2 not allowed
@@ -1059,7 +1060,7 @@ def test_fill_with_sequence_2():
         bh.axis.StrCategory(["A", "B"]),
     )
     b.fill((1, 0, 10), ("C", "B", "A"))
-    assert_array_equal(b.view(True), [[0, 1, 0], [0, 0, 1]])
+    assert np.asarray(b.view(True)) == approx(np.asarray([[0, 1, 0], [0, 0, 1]]))
 
 
 def test_fill_with_sequence_3():
@@ -1070,7 +1071,7 @@ def test_fill_with_sequence_3():
     assert h.axes[0].size == 1
     h.fill(["A", "B"])
     assert h.axes[0].size == 2
-    assert_array_equal(h.view(True), [3, 1])
+    assert np.asarray(h.view(True)) == approx(np.asarray([3, 1]))
 
 
 @pytest.mark.skipif(
@@ -1084,7 +1085,7 @@ def test_fill_with_sequence_4():
     h.fill("1", np.arange(2))
     assert h.axes[0].size == 1
     assert h.axes[1].size == 2
-    assert_array_equal(h.view(True), [[1, 1]])
+    assert np.asarray(h.view(True)) == approx(np.asarray([[1, 1]]))
 
     with pytest.raises(ValueError):
         h.fill(["1"], np.arange(2))  # lengths do not match
@@ -1168,7 +1169,7 @@ def test_hist_division():
 
     h1 /= h.axes[0].widths * h.sum()
 
-    assert_array_equal(h1.view(), dens)
+    assert np.asarray(h1.view()) == approx(np.asarray(dens))
 
 
 # issue #416 b
@@ -1184,7 +1185,7 @@ def test_hist_division():
 #
 #    h1[:] /=  h.axes[0].widths * h.sum()
 #
-#    assert_allclose(h1.view(), dens)
+#    assert np.asarray(h1.view()) == approx(np.asarray(dens)) #not sure but this should be updated in the issue too
 
 
 def test_add_hists():
@@ -1203,10 +1204,10 @@ def test_add_hists():
     h3 = h.copy()
     h3 += 5
 
-    assert_array_equal(h, 1)
-    assert_array_equal(h1, 2)
-    assert_array_equal(h2, 3)
-    assert_array_equal(h3, 6)
+    assert h == approx(1)
+    assert h1 == approx(2)
+    assert h2 == approx(3)
+    assert h3 == approx(6)
 
 
 def test_add_broadcast():
@@ -1258,7 +1259,7 @@ def test_reductions():
     widths_1 = functools.reduce(operator.mul, h.axes.widths)
     widths_2 = np.prod(h.axes.widths, axis=0)
 
-    assert_array_equal(widths_1, widths_2)
+    assert widths_1 == approx(widths_2)
 
 
 # Issue 435
diff --git a/tests/test_histogram_indexing.py b/tests/test_histogram_indexing.py
index 8206356f3..c19e76863 100644
--- a/tests/test_histogram_indexing.py
+++ b/tests/test_histogram_indexing.py
@@ -1,6 +1,5 @@
 import numpy as np
 import pytest
-from numpy.testing import assert_array_equal
 from pytest import approx
 
 import boost_histogram as bh
@@ -181,8 +180,8 @@ def test_mix_value_with_slice():
     assert h[1, 1, True] == 11
     assert h[3, 4, False] == 0
 
-    assert_array_equal(h[:, :, True].view(), vals[:, :, 0])
-    assert_array_equal(h[:, :, False].view(), 0)
+    assert np.asarray(h[:, :, True].view()) == approx(np.asarray(vals[:, :, 0]))
+    assert np.asarray(h[:, :, False].view()) == approx(np.asarray(0))
 
 
 def test_mix_value_with_slice_2():
@@ -198,11 +197,11 @@ def test_mix_value_with_slice_2():
     assert h[1, 1, True] == 11
     assert h[3, 4, False] == 0
 
-    assert_array_equal(h[:, :, True].view(), vals)
-    assert_array_equal(h[:, :, False].view(), 0)
+    assert np.asarray(h[:, :, True].view()) == approx(np.asarray(vals))
+    assert np.asarray(h[:, :, False].view()) == approx(np.asarray(0))
 
     h2 = h[bh.rebin(2), bh.rebin(5), :]
-    assert_array_equal(h2.shape, (5, 2, 2))
+    assert h2.shape == approx((5, 2, 2))
 
 
 def test_one_sided_slice():
@@ -213,7 +212,7 @@ def test_one_sided_slice():
     assert h[bh.tag.at(-1) : bh.tag.at(5) : sum] == 6  # keeps underflow, keeps overflow
 
     # check that slicing without bh.sum adds removed counts to flow bins
-    assert_array_equal(h[1:3].view(True), [2, 1, 1, 2])
+    assert np.asarray(h[1:3].view(True)) == approx(np.asarray([2, 1, 1, 2]))
 
     assert h[0::sum] == 5  # removes underflow, keeps overflow
     assert h[:4:sum] == 5  # removes overflow, keeps underflow
@@ -263,8 +262,8 @@ def test_noflow_slicing():
     assert h[3, 4, False] == 0
     assert h[{0: 3, 1: 4, 2: False}] == 0
 
-    assert_array_equal(h[:, :, True].view(), vals)
-    assert_array_equal(h[:, :, False].view(), 0)
+    assert np.asarray(h[:, :, True].view()) == approx(np.asarray(vals))
+    assert np.asarray(h[:, :, False].view()) == approx(np.asarray(0))
 
 
 def test_singleflow_slicing():
@@ -280,9 +279,9 @@ def test_singleflow_slicing():
     assert h[1, 0] == 4
     assert h[1, 1] == 5
 
-    assert_array_equal(h[:, 1 : 3 : bh.sum], vals[:, 1:3].sum(axis=1))
-    assert_array_equal(h[{1: slice(1, 3, bh.sum)}], vals[:, 1:3].sum(axis=1))
-    assert_array_equal(h[1 : 3 : bh.sum, :], vals[1:3, :].sum(axis=0))
+    assert h[:, 1 : 3 : bh.sum] == approx(vals[:, 1:3].sum(axis=1))
+    assert h[{1: slice(1, 3, bh.sum)}] == approx(vals[:, 1:3].sum(axis=1))
+    assert h[1 : 3 : bh.sum, :] == approx(vals[1:3, :].sum(axis=0))
 
 
 def test_pick_str_category():
@@ -302,9 +301,9 @@ def test_pick_str_category():
     assert h[1, 1, bh.loc("on")] == 11
     assert h[3, 4, bh.loc("maybe")] == 0
 
-    assert_array_equal(h[:, :, bh.loc("on")].view(), vals)
-    assert_array_equal(h[{2: bh.loc("on")}].view(), vals)
-    assert_array_equal(h[:, :, bh.loc("off")].view(), 0)
+    assert np.asarray(h[:, :, bh.loc("on")].view()) == approx(np.asarray(vals))
+    assert np.asarray(h[{2: bh.loc("on")}].view()) == approx(np.asarray(vals))
+    assert np.asarray(h[:, :, bh.loc("off")].view()) == approx(np.asarray(0))
 
 
 def test_string_requirement():
@@ -346,10 +345,10 @@ def test_pick_int_category():
     assert h[3, 4, bh.loc(7)] == 0
     assert h[3, 4, bh.loc(12)] == 134
 
-    assert_array_equal(h[:, :, bh.loc(3)].view(), vals)
-    assert_array_equal(h[{2: bh.loc(3)}].view(), vals)
-    assert_array_equal(h[:, :, bh.loc(5)].view(), vals + 1)
-    assert_array_equal(h[:, :, bh.loc(7)].view(), 0)
+    assert np.asarray(h[:, :, bh.loc(3)].view()) == approx(np.asarray(vals))
+    assert np.asarray(h[{2: bh.loc(3)}].view()) == approx(np.asarray(vals))
+    assert np.asarray(h[:, :, bh.loc(5)].view()) == approx(np.asarray(vals + 1))
+    assert np.asarray(h[:, :, bh.loc(7)].view()) == approx(np.asarray(0))
 
 
 @pytest.mark.parametrize(
@@ -384,7 +383,7 @@ def test_axes_tuple():
     (before,) = h.axes.centers[:1]
     (after,) = h.axes[:1].centers
 
-    assert_array_equal(before, after)
+    assert before == approx(after)
 
 
 def test_axes_tuple_Nd():
@@ -397,8 +396,8 @@ def test_axes_tuple_Nd():
     b1, b2 = h.axes.centers[1:3]
     a1, a2 = h.axes[1:3].centers
 
-    assert_array_equal(b1.flatten(), a1.flatten())
-    assert_array_equal(b2.flatten(), a2.flatten())
+    assert b1.flatten() == approx(a1.flatten())
+    assert b2.flatten() == approx(a2.flatten())
 
     assert b1.ndim == 3
     assert a1.ndim == 2
diff --git a/tests/test_histogram_set.py b/tests/test_histogram_set.py
index 92221d0c5..fbf335d4e 100644
--- a/tests/test_histogram_set.py
+++ b/tests/test_histogram_set.py
@@ -1,6 +1,6 @@
 import numpy as np
 import pytest
-from numpy.testing import assert_array_equal
+from pytest import approx
 
 import boost_histogram as bh
 
@@ -37,10 +37,10 @@ def test_1d_set_array():
     h = bh.Histogram(bh.axis.Regular(10, 0, 1))
 
     h[...] = np.arange(10)
-    assert_array_equal(h.view(), np.arange(10))
+    assert np.asarray(h.view()) == approx(np.asarray(np.arange(10)))
 
     h[...] = np.arange(12)
-    assert_array_equal(h.view(flow=True), np.arange(12))
+    assert np.asarray(h.view(flow=True)) == approx(np.asarray(np.arange(12)))
 
     with pytest.raises(ValueError):
         h[...] = np.arange(9)
@@ -50,17 +50,17 @@ def test_1d_set_array():
         h[...] = np.arange(13)
 
     h[...] = 1
-    assert_array_equal(h.view(), np.ones(10))
+    assert np.asarray(h.view()) == approx(np.asarray(np.ones(10)))
 
 
 def test_2d_set_array():
     h = bh.Histogram(bh.axis.Regular(10, 0, 1), bh.axis.Regular(10, 0, 1))
 
     h[...] = np.arange(10).reshape(-1, 1)
-    assert_array_equal(h.view()[:, 2], np.arange(10))
+    assert np.asarray(h.view()[:, 2]) == approx(np.asarray(np.arange(10)))
 
     h[...] = np.arange(12).reshape(-1, 1)
-    assert_array_equal(h.view(flow=True)[:, 3], np.arange(12))
+    assert np.asarray(h.view(flow=True)[:, 3]) == approx(np.asarray(np.arange(12)))
 
     with pytest.raises(ValueError):
         h[...] = np.arange(9).reshape(-1, 1)
@@ -70,7 +70,7 @@ def test_2d_set_array():
         h[...] = np.arange(13).reshape(-1, 1)
 
     h[...] = 1
-    assert_array_equal(h.view(), np.ones((10, 10)))
+    assert np.asarray(h.view()) == approx(np.asarray(np.ones((10, 10))))
 
 
 def test_weighted_set_shortcut():
@@ -106,27 +106,27 @@ def test_set_special_dtype(storage, default):
 
     arr = np.full((10, 1), default)
     h[...] = arr
-    assert_array_equal(h.view()[:, 1:2], arr)
+    assert np.asarray(h.view()[:, 1:2]) == approx(np.asarray(arr))
 
     arr = np.full((12, 1), default)
     h[...] = arr
-    assert_array_equal(h.view(flow=True)[:, 2:3], arr)
+    assert np.asarray(h.view(flow=True)[:, 2:3]) == approx(np.asarray(arr))
 
     arr = np.full((10, 10), default)
     h[...] = arr
-    assert_array_equal(h.view(), arr)
+    assert np.asarray(h.view()) == approx(np.asarray(arr))
 
     arr = np.full((10, 12), default)
     h[...] = arr
-    assert_array_equal(h.view(flow=True)[1:11, :], arr)
+    assert np.asarray(h.view(flow=True)[1:11, :]) == approx(np.asarray(arr))
 
     arr = np.full((12, 10), default)
     h[...] = arr
-    assert_array_equal(h.view(flow=True)[:, 1:11], arr)
+    assert np.asarray(h.view(flow=True)[:, 1:11]) == approx(np.asarray(arr))
 
     arr = np.full((12, 12), default)
     h[...] = arr
-    assert_array_equal(h.view(flow=True), arr)
+    assert np.asarray(h.view(flow=True)) == approx(np.asarray(arr))
 
     with pytest.raises(ValueError):
         arr = np.full((9, 1), default)
diff --git a/tests/test_internal_histogram.py b/tests/test_internal_histogram.py
index 3a1489cf4..47610039e 100644
--- a/tests/test_internal_histogram.py
+++ b/tests/test_internal_histogram.py
@@ -1,6 +1,5 @@
 import numpy as np
 import pytest
-from numpy.testing import assert_allclose, assert_array_equal
 from pytest import approx
 
 import boost_histogram as bh
@@ -32,9 +31,9 @@ def test_1D_fill_int(storage):
 
     H = np.array([0, 1, 2, 0, 0, 0, 0, 0, 0, 0])
 
-    assert_array_equal(np.asarray(hist), H)
-    assert_array_equal(hist.view(flow=False), H)
-    assert_array_equal(hist.view(flow=True)[1:-1], H)
+    assert np.asarray(hist) == approx(H)
+    assert np.asarray(hist.view(flow=False)) == approx(np.asarray(H))
+    assert np.asarray(hist.view(flow=True)[1:-1]) == approx(np.asarray(H))
 
     assert hist.axes[0].size == bins
     assert hist.axes[0].extent == bins + 2
@@ -57,9 +56,9 @@ def test_2D_fill_int(storage):
 
     H = np.histogram2d(*vals, bins=bins, range=ranges)[0]
 
-    assert_array_equal(np.asarray(hist), H)
-    assert_array_equal(hist.view(flow=True)[1:-1, 1:-1], H)
-    assert_array_equal(hist.view(flow=False), H)
+    assert np.asarray(hist) == approx(H)
+    assert np.asarray(hist.view(flow=True)[1:-1, 1:-1]) == approx(np.asarray(H))
+    assert np.asarray(hist.view(flow=False)) == approx(np.asarray(H))
 
     assert hist.axes[0].size == bins[0]
     assert hist.axes[0].extent == bins[0] + 2
@@ -76,9 +75,9 @@ def test_edges_histogram():
     hist.fill(vals)
 
     bins = np.asarray(hist)
-    assert_array_equal(bins, [0, 2, 2])
-    assert_array_equal(hist.view(flow=True), [0, 0, 2, 2, 0])
-    assert_array_equal(hist.view(flow=False), [0, 2, 2])
+    assert bins == approx([0, 2, 2])
+    assert np.asarray(hist.view(flow=True)) == approx(np.asarray([0, 0, 2, 2, 0]))
+    assert np.asarray(hist.view(flow=False)) == approx(np.asarray([0, 2, 2]))
 
 
 def test_int_histogram():
@@ -88,9 +87,9 @@ def test_int_histogram():
     hist.fill(vals)
 
     bins = np.asarray(hist)
-    assert_array_equal(bins, [1, 1, 1, 1])
-    assert_array_equal(hist.view(flow=False), [1, 1, 1, 1])
-    assert_array_equal(hist.view(flow=True), [2, 1, 1, 1, 1, 3])
+    assert bins == approx([1, 1, 1, 1])
+    assert np.asarray(hist.view(flow=False)) == approx(np.asarray([1, 1, 1, 1]))
+    assert np.asarray(hist.view(flow=True)) == approx(np.asarray([2, 1, 1, 1, 1, 3]))
 
 
 def test_str_categories_histogram():
@@ -131,9 +130,9 @@ def test_numpy_dd():
     h2, x2, y2 = h.to_numpy()
     h1, (x1, y1) = h.to_numpy(dd=True)
 
-    assert_array_equal(h1, h2)
-    assert_array_equal(x1, x2)
-    assert_array_equal(y1, y2)
+    assert h1 == approx(h2)
+    assert x1 == approx(x2)
+    assert y1 == approx(y2)
 
 
 def test_numpy_weights():
@@ -150,16 +149,16 @@ def test_numpy_weights():
     h2, x2, y2 = h.to_numpy(view=False)
     h1, (x1, y1) = h.to_numpy(dd=True, view=False)
 
-    assert_array_equal(h1, h2)
-    assert_array_equal(x1, x2)
-    assert_array_equal(y1, y2)
+    assert h1 == approx(h2)
+    assert x1 == approx(x2)
+    assert y1 == approx(y2)
 
     h1, (x1, y1) = h.to_numpy(dd=True, view=False)
     h2, x2, y2 = h.to_numpy(view=True)
 
-    assert_array_equal(h1, h2.value)
-    assert_array_equal(x1, x2)
-    assert_array_equal(y1, y2)
+    assert h1 == approx(h2.value)
+    assert x1 == approx(x2)
+    assert y1 == approx(y2)
 
 
 def test_numpy_flow():
@@ -176,14 +175,14 @@ def test_numpy_flow():
     flow_true = h.to_numpy(True)[0][1:-1, 1:-1]
     flow_false = h.to_numpy(False)[0]
 
-    assert_array_equal(flow_true, flow_false)
+    assert flow_true == approx(flow_false)
 
     view_flow_true = h.view(flow=True)
     view_flow_false = h.view(flow=False)
     view_flow_default = h.view()
 
-    assert_array_equal(view_flow_true[1:-1, 1:-1], view_flow_false)
-    assert_array_equal(view_flow_default, view_flow_false)
+    assert view_flow_true[1:-1, 1:-1] == approx(view_flow_false)
+    assert view_flow_default == approx(view_flow_false)
 
 
 def test_numpy_compare():
@@ -206,9 +205,9 @@ def test_numpy_compare():
 
     nH, nE1, nE2 = np.histogram2d(xs, ys, bins=(10, 5), range=((0, 1), (0, 1)))
 
-    assert_array_equal(H, nH)
-    assert_allclose(E1, nE1)
-    assert_allclose(E2, nE2)
+    assert H == approx(nH)
+    assert E1 == approx(nE1)
+    assert E2 == approx(nE2)
 
 
 def test_project():
@@ -234,9 +233,9 @@ def test_project():
     assert h.project(0) == h0
     assert h.project(1) == h1
 
-    assert_array_equal(h.project(0, 1), h)
-    assert_array_equal(h.project(0), h0)
-    assert_array_equal(h.project(1), h1)
+    assert h.project(0, 1) == approx(h)
+    assert h.project(0) == approx(h0)
+    assert h.project(1) == approx(h1)
 
 
 def test_sums():
@@ -254,7 +253,7 @@ def test_int_cat_hist():
     h.fill(2)
     h.fill(3)
 
-    assert_array_equal(h.view(), [1, 1, 1])
+    assert np.asarray(h.view()) == approx(np.asarray([1, 1, 1]))
     assert h.sum() == 3
 
     with pytest.raises(RuntimeError):
diff --git a/tests/test_plottable_protocol.py b/tests/test_plottable_protocol.py
index 5044fb488..e2933c8d2 100644
--- a/tests/test_plottable_protocol.py
+++ b/tests/test_plottable_protocol.py
@@ -1,5 +1,5 @@
 import numpy as np
-from numpy.testing import assert_allclose
+from pytest import approx
 
 import boost_histogram as bh
 
@@ -14,9 +14,9 @@ def test_plottable_histogram_mean_int():
 
     h[...] = np.stack([COUNTS, VALUES, VARIANCES]).T
 
-    assert_allclose(COUNTS, h.counts())
-    assert_allclose(VALUES, h.values())
-    assert_allclose(VARIANCES / COUNTS, h.variances())
+    assert COUNTS == approx(h.counts())
+    assert VALUES == approx(h.values())
+    assert (VARIANCES / COUNTS) == approx(h.variances())
 
     assert h.kind == bh.Kind.MEAN
     assert h.kind == "MEAN"
@@ -38,9 +38,9 @@ def test_plottible_histogram_weight_reg():
 
     h[...] = np.stack([VALUES, VARIANCES]).T
 
-    assert_allclose(VALUES, h.counts())
-    assert_allclose(VALUES, h.values())
-    assert_allclose(VARIANCES, h.variances())
+    assert VALUES == approx(h.counts())
+    assert VALUES == approx(h.values())
+    assert VARIANCES == approx(h.variances())
 
     assert h.kind == bh.Kind.COUNT
     assert h.kind == "COUNT"
@@ -48,10 +48,10 @@ def test_plottible_histogram_weight_reg():
     assert len(h.axes) == 1
     assert h.axes[0] == h.axes[0]
 
-    assert_allclose(h.axes[0][0], (0, 1))
-    assert_allclose(h.axes[0][1], (1, 2))
-    assert_allclose(h.axes[0][2], (2, 3))
-    assert_allclose(h.axes[0][3], (3, 4))
+    assert h.axes[0][0] == approx((0, 1))
+    assert h.axes[0][1] == approx((1, 2))
+    assert h.axes[0][2] == approx((2, 3))
+    assert h.axes[0][3] == approx((3, 4))
 
     assert not h.axes[0].traits.discrete
     assert not h.axes[0].traits.circular
@@ -66,9 +66,9 @@ def test_plottible_histogram_simple_var():
     # are doing. At least if you change it inplace.
     h.view()[...] = VALUES
 
-    assert_allclose(VALUES, h.counts())
-    assert_allclose(VALUES, h.values())
-    assert_allclose(VALUES, h.variances())
+    assert VALUES == approx(h.counts())
+    assert VALUES == approx(h.values())
+    assert VALUES == approx(h.variances())
 
     assert h.kind == bh.Kind.COUNT
     assert h.kind == "COUNT"
@@ -76,18 +76,18 @@ def test_plottible_histogram_simple_var():
     assert len(h.axes) == 1
     assert h.axes[0] == h.axes[0]
 
-    assert_allclose(h.axes[0][0], (0, 1))
-    assert_allclose(h.axes[0][1], (1, 2))
-    assert_allclose(h.axes[0][2], (2, 3))
-    assert_allclose(h.axes[0][3], (3, 4))
+    assert h.axes[0][0] == approx((0, 1))
+    assert h.axes[0][1] == approx((1, 2))
+    assert h.axes[0][2] == approx((2, 3))
+    assert h.axes[0][3] == approx((3, 4))
 
     assert not h.axes[0].traits.discrete
     assert not h.axes[0].traits.circular
 
     h.fill([1], weight=0)
 
-    assert_allclose(VALUES, h.counts())
-    assert_allclose(VALUES, h.values())
+    assert VALUES == approx(h.counts())
+    assert VALUES == approx(h.values())
     assert h.variances() is None
 
 
@@ -97,16 +97,16 @@ def test_plottible_histogram_simple_var_invalidate_inplace():
 
     h2 = h * 1
 
-    assert_allclose(VALUES, h.counts())
-    assert_allclose(VALUES, h.values())
-    assert_allclose(VALUES, h.variances())
+    assert VALUES == approx(h.counts())
+    assert VALUES == approx(h.values())
+    assert VALUES == approx(h.variances())
 
-    assert_allclose(VALUES, h2.counts())
-    assert_allclose(VALUES, h2.values())
+    assert VALUES == approx(h2.counts())
+    assert VALUES == approx(h2.values())
     assert h2.variances() is None
 
     h *= 1
 
-    assert_allclose(VALUES, h.counts())
-    assert_allclose(VALUES, h.values())
+    assert VALUES == approx(h.counts())
+    assert VALUES == approx(h.values())
     assert h.variances() is None
diff --git a/tests/test_storage.py b/tests/test_storage.py
index b151fabcc..010c47b51 100644
--- a/tests/test_storage.py
+++ b/tests/test_storage.py
@@ -2,7 +2,6 @@
 
 import numpy as np
 import pytest
-from numpy.testing import assert_array_equal
 from pytest import approx
 
 import boost_histogram as bh
@@ -23,7 +22,7 @@ def test_setting(storage):
     assert h[1] == 3
     assert h[9] == 5
 
-    assert_array_equal(h.view(), [2, 3, 0, 0, 0, 0, 0, 0, 0, 5])
+    assert np.asarray(h.view()) == approx(np.asarray([2, 3, 0, 0, 0, 0, 0, 0, 0, 5]))
 
 
 def test_setting_weight():
@@ -62,8 +61,8 @@ def test_setting_weight():
     assert b["value"][0] == a["value"][0]
     assert b["variance"][0] == a["variance"][0]
 
-    assert_array_equal(a.view().value, b.view()["value"])
-    assert_array_equal(a.view().variance, b.view()["variance"])
+    assert np.asarray(a.view().value) == approx(np.asarray(b.view()["value"]))
+    assert np.asarray(a.view().variance) == approx(np.asarray(b.view()["variance"]))
 
 
 def test_sum_weight():
@@ -78,7 +77,7 @@ def test_sum_weight():
     v2 = v + v
     h2 = h + h
 
-    assert_array_equal(h2.view(), v2)
+    assert np.asarray(h2.view()) == approx(np.asarray(v2))
 
 
 def test_setting_profile():
@@ -121,10 +120,10 @@ def test_setting_profile():
     assert b[0]["count"] == a["count"][0]
     assert b[0]["_sum_of_deltas_squared"] == a["_sum_of_deltas_squared"][0]
 
-    assert_array_equal(a.view().value, b.view()["value"])
-    assert_array_equal(a.view().count, b.view()["count"])
-    assert_array_equal(
-        a.view()._sum_of_deltas_squared, b.view()["_sum_of_deltas_squared"]
+    assert np.asarray(a.view().value) == approx(np.asarray(b.view()["value"]))
+    assert np.asarray(a.view().count) == approx(np.asarray(b.view()["count"]))
+    assert np.asarray(a.view()._sum_of_deltas_squared) == approx(
+        np.asarray(b.view()["_sum_of_deltas_squared"])
     )
 
 
@@ -195,14 +194,17 @@ def test_setting_weighted_profile():
         == a["_sum_of_weighted_deltas_squared"][0]
     )
 
-    assert_array_equal(a.view().value, b.view()["value"])
-    assert_array_equal(a.view().sum_of_weights, b.view()["sum_of_weights"])
-    assert_array_equal(
-        a.view().sum_of_weights_squared, b.view()["sum_of_weights_squared"]
+    assert np.asarray(a.view().value) == approx(np.asarray(b.view()["value"]))
+    assert np.asarray(a.view().sum_of_weights) == approx(
+        np.asarray(b.view()["sum_of_weights"])
     )
-    assert_array_equal(
-        a.view()._sum_of_weighted_deltas_squared,
-        b.view()["_sum_of_weighted_deltas_squared"],
+    assert np.asarray(a.view().sum_of_weights_squared) == approx(
+        np.asarray(b.view()["sum_of_weights_squared"])
+    )
+    assert np.asarray(a.view()._sum_of_weighted_deltas_squared) == approx(
+        np.asarray(
+            b.view()["_sum_of_weighted_deltas_squared"],
+        )
     )
 
 
@@ -216,8 +218,8 @@ def test_modify_weights_by_view():
 
     hist.view().value /= 2
 
-    assert hist.view().value[0] == pytest.approx(1.5)
-    assert hist.view().value[1] == pytest.approx(2)
+    assert hist.view().value[0] == approx(1.5)
+    assert hist.view().value[1] == approx(2)
 
 
 # Issue #531
diff --git a/tests/test_threaded_fill.py b/tests/test_threaded_fill.py
index 627916aa5..1921a59e5 100644
--- a/tests/test_threaded_fill.py
+++ b/tests/test_threaded_fill.py
@@ -1,6 +1,6 @@
 import numpy as np
 import pytest
-from numpy.testing import assert_almost_equal, assert_array_equal
+from pytest import approx
 
 import boost_histogram as bh
 
@@ -24,7 +24,7 @@ def test_threads(benchmark, threads, storage):
     hist_linear.fill(vals)
     hist_result = benchmark(fillit, hist_atomic, vals, threads=threads)
 
-    assert_array_equal(hist_linear, hist_result)
+    assert hist_linear == approx(hist_result)
 
 
 @pytest.mark.parametrize("threads", [1, 4, 7], ids=lambda x: f"threads={x}")
@@ -41,7 +41,7 @@ def test_threaded_builtin(threads, storage):
     hist_atomic1.fill(vals)
     hist_atomic2.fill(vals, threads=threads)
 
-    assert_array_equal(hist_atomic1, hist_atomic2)
+    assert hist_atomic1 == approx(hist_atomic2)
 
 
 @pytest.mark.parametrize("threads", [1, 4, 7], ids=lambda x: f"threads={x}")
@@ -51,7 +51,7 @@ def test_threaded_numpy(threads):
     hist_1, _ = bh.numpy.histogram(vals)
     hist_2, _ = bh.numpy.histogram(vals, threads=threads)
 
-    assert_array_equal(hist_1, hist_2)
+    assert hist_1 == approx(hist_2)
 
 
 @pytest.mark.parametrize("threads", [1, 4, 7], ids=lambda x: f"threads={x}")
@@ -64,7 +64,7 @@ def test_threaded_weights(threads):
     hist_1.fill(x, y, weight=weights)
     hist_2.fill(x, y, weight=weights, threads=threads)
 
-    assert_almost_equal(hist_1.view(), hist_2.view())
+    assert np.asarray(hist_1.view()) == approx(np.asarray(hist_2.view()))
 
 
 @pytest.mark.parametrize("threads", [1, 4, 7], ids=lambda x: f"threads={x}")
@@ -81,8 +81,10 @@ def test_threaded_weight_storage(threads):
     hist_1.fill(x, y, weight=weights)
     hist_2.fill(x, y, weight=weights, threads=threads)
 
-    assert_almost_equal(hist_1.view().value, hist_2.view().value)
-    assert_almost_equal(hist_1.view().variance, hist_2.view().variance)
+    assert np.asarray(hist_1.view().value) == approx(np.asarray(hist_2.view().value))
+    assert np.asarray(hist_1.view().variance) == approx(
+        np.asarray(hist_2.view().variance)
+    )
 
 
 def test_no_profile():
diff --git a/tests/test_views.py b/tests/test_views.py
index e70b63942..6464c9a27 100644
--- a/tests/test_views.py
+++ b/tests/test_views.py
@@ -1,6 +1,5 @@
 import numpy as np
 import pytest
-from numpy.testing import assert_allclose
 from pytest import approx
 
 import boost_histogram as bh
@@ -14,115 +13,115 @@ def v():
 
 
 def test_basic_view(v):
-    assert_allclose(v.value, [0, 3, 2, 1])
-    assert_allclose(v.variance, [0, 3, 2, 1])
+    assert np.asarray(v.value) == approx(np.asarray([0, 3, 2, 1]))
+    assert np.asarray(v.variance) == approx(np.asarray([0, 3, 2, 1]))
 
 
 def test_view_mul(v):
     v2 = v * 2
-    assert_allclose(v2.value, [0, 6, 4, 2])
-    assert_allclose(v2.variance, [0, 12, 8, 4])
+    assert np.asarray(v2.value) == approx(np.asarray([0, 6, 4, 2]))
+    assert np.asarray(v2.variance) == approx(np.asarray([0, 12, 8, 4]))
 
     v2 = 2 * v
-    assert_allclose(v2.value, [0, 6, 4, 2])
-    assert_allclose(v2.variance, [0, 12, 8, 4])
+    assert np.asarray(v2.value) == approx(np.asarray([0, 6, 4, 2]))
+    assert np.asarray(v2.variance) == approx(np.asarray([0, 12, 8, 4]))
 
     v2 = v * (-2)
-    assert_allclose(v2.value, [0, -6, -4, -2])
-    assert_allclose(v2.variance, [0, 12, 8, 4])
+    assert np.asarray(v2.value) == approx(np.asarray([0, -6, -4, -2]))
+    assert np.asarray(v2.variance) == approx(np.asarray([0, 12, 8, 4]))
 
     v *= 2
-    assert_allclose(v.value, [0, 6, 4, 2])
-    assert_allclose(v.variance, [0, 12, 8, 4])
+    assert np.asarray(v.value) == approx(np.asarray([0, 6, 4, 2]))
+    assert np.asarray(v.variance) == approx(np.asarray([0, 12, 8, 4]))
 
 
 def test_view_div(v):
     v2 = v / 2
-    assert_allclose(v2.value, [0, 1.5, 1, 0.5])
-    assert_allclose(v2.variance, [0, 0.75, 0.5, 0.25])
+    assert np.asarray(v2.value) == approx(np.asarray([0, 1.5, 1, 0.5]))
+    assert np.asarray(v2.variance) == approx(np.asarray([0, 0.75, 0.5, 0.25]))
 
     v2 = v / (-0.5)
-    assert_allclose(v2.value, [0, -6, -4, -2])
-    assert_allclose(v2.variance, [0, 12, 8, 4])
+    assert np.asarray(v2.value) == approx(np.asarray([0, -6, -4, -2]))
+    assert np.asarray(v2.variance) == approx(np.asarray([0, 12, 8, 4]))
 
     v2 = 1 / v[1:]
-    assert_allclose(v2.value, [1 / 3, 1 / 2, 1])
-    assert_allclose(v2.variance, [1 / 3, 1 / 2, 1])
+    assert np.asarray(v2.value) == approx(np.asarray([1 / 3, 1 / 2, 1]))
+    assert np.asarray(v2.variance) == approx(np.asarray([1 / 3, 1 / 2, 1]))
 
     v /= 0.5
-    assert_allclose(v.value, [0, 6, 4, 2])
-    assert_allclose(v.variance, [0, 12, 8, 4])
+    assert np.asarray(v.value) == approx(np.asarray([0, 6, 4, 2]))
+    assert np.asarray(v.variance) == approx(np.asarray([0, 12, 8, 4]))
 
 
 def test_view_add(v):
     v2 = v + 1
-    assert_allclose(v2.value, [1, 4, 3, 2])
-    assert_allclose(v2.variance, [1, 4, 3, 2])
+    assert np.asarray(v2.value) == approx(np.asarray([1, 4, 3, 2]))
+    assert np.asarray(v2.variance) == approx(np.asarray([1, 4, 3, 2]))
 
     v2 = v + 2
-    assert_allclose(v2.value, [2, 5, 4, 3])
-    assert_allclose(v2.variance, [4, 7, 6, 5])
+    assert np.asarray(v2.value) == approx(np.asarray([2, 5, 4, 3]))
+    assert np.asarray(v2.variance) == approx(np.asarray([4, 7, 6, 5]))
 
     v2 = 2 + v
-    assert_allclose(v2.value, [2, 5, 4, 3])
-    assert_allclose(v2.variance, [4, 7, 6, 5])
+    assert np.asarray(v2.value) == approx(np.asarray([2, 5, 4, 3]))
+    assert np.asarray(v2.variance) == approx(np.asarray([4, 7, 6, 5]))
 
     v2 = v.copy()
     v2 += 2
-    assert_allclose(v2.value, [2, 5, 4, 3])
-    assert_allclose(v2.variance, [4, 7, 6, 5])
+    assert np.asarray(v2.value) == approx(np.asarray([2, 5, 4, 3]))
+    assert np.asarray(v2.variance) == approx(np.asarray([4, 7, 6, 5]))
 
     v2 = v + v
-    assert_allclose(v2.value, v.value * 2)
-    assert_allclose(v2.variance, v.variance * 2)
+    assert np.asarray(v2.value) == approx(np.asarray(v.value * 2))
+    assert np.asarray(v2.variance) == approx(np.asarray(v.variance * 2))
 
 
 def test_view_sub(v):
     v2 = v - 1
-    assert_allclose(v2.value, [-1, 2, 1, 0])
-    assert_allclose(v2.variance, [1, 4, 3, 2])
+    assert np.asarray(v2.value) == approx(np.asarray([-1, 2, 1, 0]))
+    assert np.asarray(v2.variance) == approx(np.asarray([1, 4, 3, 2]))
 
     v2 = v - 2
-    assert_allclose(v2.value, [-2, 1, 0, -1])
-    assert_allclose(v2.variance, [4, 7, 6, 5])
+    assert np.asarray(v2.value) == approx(np.asarray([-2, 1, 0, -1]))
+    assert np.asarray(v2.variance) == approx(np.asarray([4, 7, 6, 5]))
 
     v2 = 1 - v
-    assert_allclose(v2.value, [1, -2, -1, 0])
-    assert_allclose(v2.variance, [1, 4, 3, 2])
+    assert np.asarray(v2.value) == approx(np.asarray([1, -2, -1, 0]))
+    assert np.asarray(v2.variance) == approx(np.asarray([1, 4, 3, 2]))
 
     v2 = v.copy()
     v2 -= 2
-    assert_allclose(v2.value, [-2, 1, 0, -1])
-    assert_allclose(v2.variance, [4, 7, 6, 5])
+    assert np.asarray(v2.value) == approx(np.asarray([-2, 1, 0, -1]))
+    assert np.asarray(v2.variance) == approx(np.asarray([4, 7, 6, 5]))
 
     v2 = v - v
-    assert_allclose(v2.value, [0, 0, 0, 0])
-    assert_allclose(v2.variance, v.variance * 2)
+    assert np.asarray(v2.value) == approx(np.asarray([0, 0, 0, 0]))
+    assert np.asarray(v2.variance) == approx(np.asarray(v.variance * 2))
 
 
 def test_view_unary(v):
     v2 = +v
-    assert_allclose(v.value, v2.value)
-    assert_allclose(v.variance, v2.variance)
+    assert np.asarray(v.value) == approx(np.asarray(v2.value))
+    assert np.asarray(v.variance) == approx(np.asarray(v2.variance))
 
     v2 = -v
-    assert_allclose(-v.value, v2.value)
-    assert_allclose(v.variance, v2.variance)
+    assert np.asarray(-v.value) == approx(np.asarray(v2.value))
+    assert np.asarray(v.variance) == approx(np.asarray(v2.variance))
 
 
 def test_view_add_same(v):
     v2 = v + v
 
-    assert_allclose(v.value * 2, v2.value)
-    assert_allclose(v.variance * 2, v2.variance)
+    assert np.asarray(v.value * 2) == approx(np.asarray(v2.value))
+    assert np.asarray(v.variance * 2) == approx(np.asarray(v2.variance))
 
     v2 = v + v[1]
-    assert_allclose(v.value + 3, v2.value)
-    assert_allclose(v.variance + 3, v2.variance)
+    assert np.asarray(v.value + 3) == approx(np.asarray(v2.value))
+    assert np.asarray(v.variance + 3) == approx(np.asarray(v2.variance))
 
     v2 = v + bh.accumulators.WeightedSum(5, 6)
-    assert_allclose(v.value + 5, v2.value)
-    assert_allclose(v.variance + 6, v2.variance)
+    assert np.asarray(v.value + 5) == approx(np.asarray(v2.value))
+    assert np.asarray(v.variance + 6) == approx(np.asarray(v2.variance))
 
     with pytest.raises(TypeError):
         v2 = v + bh.accumulators.WeightedMean(1, 2, 5, 6)
@@ -131,8 +130,8 @@ def test_view_add_same(v):
 def test_view_assign(v):
     v[...] = [[4, 1], [5, 2], [6, 1], [7, 2]]
 
-    assert_allclose(v.value, [4, 5, 6, 7])
-    assert_allclose(v.variance, [1, 2, 1, 2])
+    assert np.asarray(v.value) == approx(np.asarray([4, 5, 6, 7]))
+    assert np.asarray(v.variance) == approx(np.asarray([1, 2, 1, 2]))
 
 
 def test_view_assign_mean():
@@ -140,9 +139,9 @@ def test_view_assign_mean():
     m = h.copy().view()
 
     h[...] = [[10, 2, 3], [4, 5, 6], [7, 8, 9], [10, 11, 12]]
-    assert_allclose(h.view().count, [10, 4, 7, 10])
-    assert_allclose(h.view().value, [2, 5, 8, 11])
-    assert_allclose(h.view().variance, [3, 6, 9, 12])
+    assert np.asarray(h.view().count) == approx(np.asarray([10, 4, 7, 10]))
+    assert np.asarray(h.view().value) == approx(np.asarray([2, 5, 8, 11]))
+    assert np.asarray(h.view().variance) == approx(np.asarray([3, 6, 9, 12]))
 
     # Make sure this really was a copy
     assert m.count[0] != 10
@@ -150,9 +149,9 @@ def test_view_assign_mean():
     # Assign directly on view
     m[...] = [[10, 2, 3], [4, 5, 6], [7, 8, 9], [10, 11, 12]]
 
-    assert_allclose(m.count, [10, 4, 7, 10])
-    assert_allclose(m.value, [2, 5, 8, 11])
-    assert_allclose(m.variance, [3, 6, 9, 12])
+    assert np.asarray(m.count) == approx(np.asarray([10, 4, 7, 10]))
+    assert np.asarray(m.value) == approx(np.asarray([2, 5, 8, 11]))
+    assert np.asarray(m.variance) == approx(np.asarray([3, 6, 9, 12]))
     # Note: if counts <= 1, variance is undefined
 
 
@@ -163,10 +162,12 @@ def test_view_assign_wmean():
 
     h[...] = [[10, 2, 3, 4], [5, 6, 7, 8], [9, 10, 11, 12], [13, 14, 15, 16]]
 
-    assert_allclose(h.view().sum_of_weights, [10, 5, 9, 13])
-    assert_allclose(h.view().sum_of_weights_squared, [2, 6, 10, 14])
-    assert_allclose(h.view().value, [3, 7, 11, 15])
-    assert_allclose(h.view().variance, [4, 8, 12, 16])
+    assert np.asarray(h.view().sum_of_weights) == approx(np.asarray([10, 5, 9, 13]))
+    assert np.asarray(h.view().sum_of_weights_squared) == approx(
+        np.asarray([2, 6, 10, 14])
+    )
+    assert np.asarray(h.view().value) == approx(np.asarray([3, 7, 11, 15]))
+    assert np.asarray(h.view().variance) == approx(np.asarray([4, 8, 12, 16]))
 
     # Make sure this really was a copy
     assert w.sum_of_weights[0] != 10
@@ -174,10 +175,10 @@ def test_view_assign_wmean():
     # Assign directly on view
     w[...] = [[10, 2, 3, 4], [5, 6, 7, 8], [9, 10, 11, 12], [13, 14, 15, 16]]
 
-    assert_allclose(w.sum_of_weights, [10, 5, 9, 13])
-    assert_allclose(w.sum_of_weights_squared, [2, 6, 10, 14])
-    assert_allclose(w.value, [3, 7, 11, 15])
-    assert_allclose(w.variance, [4, 8, 12, 16])
+    assert np.asarray(w.sum_of_weights) == approx(np.asarray([10, 5, 9, 13]))
+    assert np.asarray(w.sum_of_weights_squared) == approx(np.asarray([2, 6, 10, 14]))
+    assert np.asarray(w.value) == approx(np.asarray([3, 7, 11, 15]))
+    assert np.asarray(w.variance) == approx(np.asarray([4, 8, 12, 16]))
     # Note: if sum_of_weights <= 1, variance is undefined
 
     w[0] = [9, 1, 2, 3]