diff --git a/doc/examples/filters/plot_inpaint.py b/doc/examples/filters/plot_inpaint.py
index ce661843..104694ec 100644
--- a/doc/examples/filters/plot_inpaint.py
+++ b/doc/examples/filters/plot_inpaint.py
@@ -44,15 +44,15 @@ image_result = inpaint.inpaint_biharmonic(image_defect, mask, multichannel=True)
 fig, axes = plt.subplots(ncols=3, nrows=1)
 
 axes[0].set_title('Defected image')
-axes[0].imshow(image_orig, cmap=plt.cm.gray, interpolation='nearest')
+axes[0].imshow(image_orig)
 axes[0].set_xticks([]), axes[0].set_yticks([])
 
 axes[1].set_title('Defect mask')
-axes[1].imshow(mask, cmap=plt.cm.gray, interpolation='nearest')
+axes[1].imshow(mask, cmap=plt.cm.gray)
 axes[1].set_xticks([]), axes[1].set_yticks([])
 
 axes[2].set_title('Inpainted image')
-axes[2].imshow(image_result, cmap=plt.cm.gray, interpolation='nearest')
+axes[2].imshow(image_result)
 axes[2].set_xticks([]), axes[2].set_yticks([])
 
 plt.show()
diff --git a/skimage/restoration/tests/test_inpaint.py b/skimage/restoration/tests/test_inpaint.py
index abe89df3..04e8c8c2 100644
--- a/skimage/restoration/tests/test_inpaint.py
+++ b/skimage/restoration/tests/test_inpaint.py
@@ -6,12 +6,13 @@ from numpy.testing import (run_module_suite, assert_allclose,
 from skimage.restoration import inpaint
 
 
-def test_inpaint_biharmonic():
+def test_inpaint_biharmonic_2d():
     img = np.tile(np.square(np.linspace(0, 1, 5)), (5, 1))
     mask = np.zeros_like(img)
     mask[2, 2:] = 1
     mask[1, 3:] = 1
     mask[0, 4:] = 1
+    img[np.where(mask)] = 0
     out = inpaint.inpaint_biharmonic(img, mask)
     ref = np.array(
         [[0., 0.0625, 0.25000000, 0.5625000, 0.73925058],
@@ -23,6 +24,32 @@ def test_inpaint_biharmonic():
     assert_allclose(ref, out)
 
 
+def test_inpaint_biharmonic_3d():
+    img = np.tile(np.square(np.linspace(0, 1, 5)), (5, 1))
+    img = np.dstack((img, img.T))
+    mask = np.zeros_like(img)
+    mask[2, 2:, :] = 1
+    mask[1, 3:, :] = 1
+    mask[0, 4:, :] = 1
+    img[np.where(mask)] = 0
+    out = inpaint.inpaint_biharmonic(img, mask)
+    ref = np.dstack((
+        np.array(
+            [[0.0000, 0.0625, 0.25000000, 0.56250000, 0.53752796],
+             [0.0000, 0.0625, 0.25000000, 0.44443780, 0.53762210],
+             [0.0000, 0.0625, 0.23693666, 0.46621112, 0.68615592],
+             [0.0000, 0.0625, 0.25000000, 0.56250000, 1.00000000],
+             [0.0000, 0.0625, 0.25000000, 0.56250000, 1.00000000]]),
+        np.array(
+            [[0.0000, 0.0000, 0.00000000, 0.00000000, 0.19621902],
+             [0.0625, 0.0625, 0.06250000, 0.17470756, 0.30140091],
+             [0.2500, 0.2500, 0.27241289, 0.35155440, 0.43068654],
+             [0.5625, 0.5625, 0.56250000, 0.56250000, 0.56250000],
+             [1.0000, 1.0000, 1.00000000, 1.00000000, 1.00000000]])
+    ))
+    assert_allclose(ref, out)
+
+
 def test_invalid_input():
     img, mask = np.zeros([]), np.zeros([])
     assert_raises(ValueError, inpaint.inpaint_biharmonic, img, mask)