diff --git a/skimage/graph/tests/test_mcp.py b/skimage/graph/tests/test_mcp.py
index a1021380..560f19d0 100644
--- a/skimage/graph/tests/test_mcp.py
+++ b/skimage/graph/tests/test_mcp.py
@@ -152,4 +152,4 @@ def _test_random(shape):
 
 
 if __name__ == "__main__":
-    run_module_suite()
+    np.testing.run_module_suite()
diff --git a/skimage/graph/tests/test_spath.py b/skimage/graph/tests/test_spath.py
index 62f9f303..d018449a 100644
--- a/skimage/graph/tests/test_spath.py
+++ b/skimage/graph/tests/test_spath.py
@@ -1,5 +1,5 @@
 import numpy as np
-from numpy.testing import *
+from numpy.testing import assert_equal, assert_array_equal
 
 import skimage.graph.spath as spath
 
@@ -33,4 +33,4 @@ def test_non_square():
 
 
 if __name__ == "__main__":
-    run_module_suite()
+    np.testing.run_module_suite()
diff --git a/skimage/measure/_find_contours.pyx b/skimage/measure/_find_contours.pyx
index d05d9aa7..e2f5a49b 100644
--- a/skimage/measure/_find_contours.pyx
+++ b/skimage/measure/_find_contours.pyx
@@ -4,8 +4,6 @@
 #cython: wraparound=False
 import numpy as np
 
-cimport numpy as cnp
-
 
 cdef inline double _get_fraction(double from_value, double to_value,
                                  double level):
@@ -14,7 +12,7 @@ cdef inline double _get_fraction(double from_value, double to_value,
     return ((level - from_value) / (to_value - from_value))
 
 
-def iterate_and_store(cnp.ndarray[double, ndim=2] array,
+def iterate_and_store(double[:, :] array,
                       double level, Py_ssize_t vertex_connect_high):
     """Iterate across the given array in a marching-squares fashion,
     looking for segments that cross 'level'. If such a segment is
@@ -46,7 +44,7 @@ def iterate_and_store(cnp.ndarray[double, ndim=2] array,
 
     # Calculate the number of iterations we'll need
     cdef Py_ssize_t num_square_steps = (array.shape[0] - 1) \
-                                       * (array.shape[1] - 1)
+                                        * (array.shape[1] - 1)
 
     cdef unsigned char square_case = 0
     cdef tuple top, bottom, left, right
diff --git a/skimage/morphology/_convex_hull.pyx b/skimage/morphology/_convex_hull.pyx
index 7298e1ed..cd9270cc 100644
--- a/skimage/morphology/_convex_hull.pyx
+++ b/skimage/morphology/_convex_hull.pyx
@@ -7,7 +7,7 @@ import numpy as np
 cimport numpy as cnp
 
 
-def possible_hull(cnp.ndarray[dtype=cnp.uint8_t, ndim=2, mode="c"] img):
+def possible_hull(cnp.uint8_t[:, ::1] img):
     """Return positions of pixels that possibly belong to the convex hull.
 
     Parameters
@@ -30,31 +30,43 @@ def possible_hull(cnp.ndarray[dtype=cnp.uint8_t, ndim=2, mode="c"] img):
     #         cols storage slots for top boundary pixels
     #         rows storage slots for right boundary pixels
     #         cols storage slots for bottom boundary pixels
-    cdef cnp.ndarray[dtype=cnp.intp_t, ndim=2] nonzero = \
-         np.ones((2 * (rows + cols), 2), dtype=np.intp)
-    nonzero *= -1
+    coords = np.ones((2 * (rows + cols), 2), dtype=np.intp)
+    coords *= -1
+
+    cdef Py_ssize_t[:, ::1] nonzero = coords
+    cdef Py_ssize_t rows_cols = rows + cols
+    cdef Py_ssize_t rows_2_cols = 2 * rows + cols
+    cdef Py_ssize_t rows_cols_r, rows_c
 
     for r in range(rows):
+
+        rows_cols_r = rows_cols + r
+
         for c in range(cols):
+
             if img[r, c] != 0:
+
+                rows_c = rows + c
+                rows_2_cols_c = rows_2_cols + c
+
                 # Left check
                 if nonzero[r, 1] == -1:
                     nonzero[r, 0] = r
                     nonzero[r, 1] = c
 
                 # Right check
-                elif nonzero[rows + cols + r, 1] < c:
-                    nonzero[rows + cols + r, 0] = r
-                    nonzero[rows + cols + r, 1] = c
+                elif nonzero[rows_cols_r, 1] < c:
+                    nonzero[rows_cols_r, 0] = r
+                    nonzero[rows_cols_r, 1] = c
 
                 # Top check
-                if nonzero[rows + c, 1] == -1:
-                    nonzero[rows + c, 0] = r
-                    nonzero[rows + c, 1] = c
+                if nonzero[rows_c, 1] == -1:
+                    nonzero[rows_c, 0] = r
+                    nonzero[rows_c, 1] = c
 
                 # Bottom check
-                elif nonzero[2 * rows + cols + c, 0] < r:
-                    nonzero[2 * rows + cols + c, 0] = r
-                    nonzero[2 * rows + cols + c, 1] = c
+                elif nonzero[rows_2_cols_c, 0] < r:
+                    nonzero[rows_2_cols_c, 0] = r
+                    nonzero[rows_2_cols_c, 1] = c
 
-    return nonzero[nonzero[:, 0] != -1]
+    return coords[coords[:, 0] != -1]
diff --git a/skimage/morphology/_greyreconstruct.pyx b/skimage/morphology/_greyreconstruct.pyx
index e8a84f3b..fa92ecec 100644
--- a/skimage/morphology/_greyreconstruct.pyx
+++ b/skimage/morphology/_greyreconstruct.pyx
@@ -21,8 +21,8 @@ def reconstruction_loop(cnp.ndarray[dtype=cnp.uint32_t, ndim=1,
                                     negative_indices=False, mode='c'] anext,
                         cnp.ndarray[dtype=cnp.int32_t, ndim=1,
                                     negative_indices=False, mode='c'] astrides,
-                        int current_idx,
-                        int image_stride):
+                        Py_ssize_t current_idx,
+                        Py_ssize_t image_stride):
     """The inner loop for reconstruction.
 
     This algorithm uses the rank-order of pixels. If low intensity pixels have
diff --git a/skimage/morphology/_pnpoly.pyx b/skimage/morphology/_pnpoly.pyx
index f32778cc..12b48e5d 100644
--- a/skimage/morphology/_pnpoly.pyx
+++ b/skimage/morphology/_pnpoly.pyx
@@ -29,7 +29,7 @@ def grid_points_inside_poly(shape, verts):
         True where the grid falls inside the polygon.
 
     """
-    cdef cnp.ndarray[cnp.double_t, ndim=1, mode="c"] vx, vy
+    cdef double[:] vx, vy
     verts = np.asarray(verts)
 
     vx = verts[:, 0].astype(np.double)
@@ -45,8 +45,7 @@ def grid_points_inside_poly(shape, verts):
 
     for m in range(M):
         for n in range(N):
-            out[m, n] = point_in_polygon(V, <double*>vx.data, <double*>vy.data,
-                                         m, n)
+            out[m, n] = point_in_polygon(V, &vx[0], &vy[0], m, n)
 
     return out.view(bool)
 
@@ -57,18 +56,18 @@ def points_inside_poly(points, verts):
     Parameters
     ----------
     points : (N, 2) array
-     Input points, ``(x, y)``.
+        Input points, ``(x, y)``.
     verts : (M, 2) array
-     Vertices of the polygon, sorted either clockwise or anti-clockwise.
-     The first point may (but does not need to be) duplicated.
+        Vertices of the polygon, sorted either clockwise or anti-clockwise.
+        The first point may (but does not need to be) duplicated.
 
     Returns
     -------
     mask : (N,) array of bool
-     True if corresponding point is inside the polygon.
+        True if corresponding point is inside the polygon.
 
     """
-    cdef cnp.ndarray[cnp.double_t, ndim=1, mode="c"] x, y, vx, vy
+    cdef double[:] x, y, vx, vy
 
     points = np.asarray(points)
     verts = np.asarray(verts)
@@ -82,8 +81,8 @@ def points_inside_poly(points, verts):
     cdef cnp.ndarray[cnp.uint8_t, ndim=1] out = \
          np.zeros(x.shape[0], dtype=np.uint8)
 
-    points_in_polygon(vx.shape[0], <double*>vx.data, <double*>vy.data,
-                      x.shape[0], <double*>x.data, <double*>y.data,
+    points_in_polygon(vx.shape[0], &vx[0], &vy[0],
+                      x.shape[0], &x[0], &y[0],
                       <unsigned char*>out.data)
 
     return out.astype(bool)
diff --git a/skimage/morphology/_skeletonize_cy.pyx b/skimage/morphology/_skeletonize_cy.pyx
index 13e303d4..10e70d98 100644
--- a/skimage/morphology/_skeletonize_cy.pyx
+++ b/skimage/morphology/_skeletonize_cy.pyx
@@ -19,16 +19,9 @@ import numpy as np
 cimport numpy as cnp
 
 
-def _skeletonize_loop(cnp.ndarray[dtype=cnp.uint8_t, ndim=2,
-                                  negative_indices=False, mode='c'] result,
-                      cnp.ndarray[dtype=cnp.intp_t, ndim=1,
-                                  negative_indices=False, mode='c'] i,
-                      cnp.ndarray[dtype=cnp.intp_t, ndim=1,
-                                  negative_indices=False, mode='c'] j,
-                      cnp.ndarray[dtype=cnp.int32_t, ndim=1,
-                                  negative_indices=False, mode='c'] order,
-                      cnp.ndarray[dtype=cnp.uint8_t, ndim=1,
-                                  negative_indices=False, mode='c'] table):
+def _skeletonize_loop(cnp.uint8_t[:, ::1] result,
+                      Py_ssize_t[:] i, Py_ssize_t[:] j,
+                      cnp.int32_t[:] order, cnp.uint8_t[:] table):
     """
     Inner loop of skeletonize function
 
@@ -65,9 +58,11 @@ def _skeletonize_loop(cnp.ndarray[dtype=cnp.uint8_t, ndim=2,
     pixels.
     """
     cdef:
-        cnp.int32_t accumulator
+        Py_ssize_t accumulator
         Py_ssize_t index, order_index
         Py_ssize_t ii, jj
+        Py_ssize_t rows = result.shape[0]
+        Py_ssize_t cols = result.shape[1]
 
     for index in range(order.shape[0]):
         accumulator = 16
@@ -80,26 +75,25 @@ def _skeletonize_loop(cnp.ndarray[dtype=cnp.uint8_t, ndim=2,
                 accumulator += 1
             if result[ii - 1, jj]:
                 accumulator += 2
-            if jj < result.shape[1] - 1 and result[ii - 1, jj + 1]:
+            if jj < cols - 1 and result[ii - 1, jj + 1]:
                     accumulator += 4
             if jj > 0 and result[ii, jj - 1]:
                 accumulator += 8
-            if jj < result.shape[1] - 1 and result[ii, jj + 1]:
+            if jj < cols - 1 and result[ii, jj + 1]:
                 accumulator += 32
-            if ii < result.shape[0]-1:
+            if ii < rows - 1:
                 if jj > 0 and result[ii + 1, jj - 1]:
                     accumulator += 64
                 if result[ii + 1, jj]:
                     accumulator += 128
-                if jj < result.shape[1] - 1 and result[ii + 1, jj + 1]:
+                if jj < cols - 1 and result[ii + 1, jj + 1]:
                     accumulator += 256
             # Assign the value of table corresponding to the configuration
             result[ii, jj] = table[accumulator]
 
 
 
-def _table_lookup_index(cnp.ndarray[dtype=cnp.uint8_t, ndim=2,
-                                    negative_indices=False, mode='c'] image):
+def _table_lookup_index(cnp.uint8_t[:, ::1] image):
     """
     Return an index into a table per pixel of a binary image
 
@@ -120,9 +114,8 @@ def _table_lookup_index(cnp.ndarray[dtype=cnp.uint8_t, ndim=2,
     hardwired kernel.
     """
     cdef:
-        cnp.ndarray[dtype=cnp.int32_t, ndim=2,
-                    negative_indices=False, mode='c'] indexer
-        cnp.int32_t *p_indexer
+        Py_ssize_t[:, ::1] indexer
+        Py_ssize_t *p_indexer
         cnp.uint8_t *p_image
         Py_ssize_t i_stride
         Py_ssize_t i_shape
@@ -133,9 +126,9 @@ def _table_lookup_index(cnp.ndarray[dtype=cnp.uint8_t, ndim=2,
 
     i_shape   = image.shape[0]
     j_shape   = image.shape[1]
-    indexer = np.zeros((i_shape, j_shape), np.int32)
-    p_indexer = <cnp.int32_t *>indexer.data
-    p_image   = <cnp.uint8_t *>image.data
+    indexer = np.zeros((i_shape, j_shape), dtype=np.intp)
+    p_indexer = &indexer[0, 0]
+    p_image   = &image[0, 0]
     i_stride  = image.strides[0]
     assert i_shape >= 3 and j_shape >= 3, \
         "Please use the slow method for arrays < 3x3"
@@ -214,4 +207,4 @@ def _table_lookup_index(cnp.ndarray[dtype=cnp.uint8_t, ndim=2,
             indexer[i - 1, j_shape - 1] += 128
             indexer[i, j_shape - 1] += 16
             indexer[i + 1, j_shape - 1] += 2
-    return indexer
+    return np.asarray(indexer)
diff --git a/skimage/morphology/_watershed.pyx b/skimage/morphology/_watershed.pyx
index 122f0262..beaf874b 100644
--- a/skimage/morphology/_watershed.pyx
+++ b/skimage/morphology/_watershed.pyx
@@ -23,17 +23,12 @@ include "heap_watershed.pxi"
 
 
 @cython.boundscheck(False)
-def watershed(np.ndarray[DTYPE_INT32_t, ndim=1, negative_indices=False,
-                         mode='c'] image,
-              np.ndarray[DTYPE_INT32_t, ndim=2, negative_indices=False,
-                         mode='c'] pq,
+def watershed(DTYPE_INT32_t[:] image,
+              DTYPE_INT32_t[:, ::1] pq,
               Py_ssize_t age,
-              np.ndarray[DTYPE_INT32_t, ndim=2, negative_indices=False,
-                         mode='c'] structure,
-              np.ndarray[DTYPE_BOOL_t, ndim=1, negative_indices=False,
-                         mode='c'] mask,
-              np.ndarray[DTYPE_INT32_t, ndim=1, negative_indices=False,
-                         mode='c'] output):
+              DTYPE_INT32_t[:, ::1] structure,
+              DTYPE_BOOL_t[:] mask,
+              DTYPE_INT32_t[:] output):
     """Do heavy lifting of watershed algorithm
 
     Parameters
diff --git a/skimage/morphology/cmorph.pyx b/skimage/morphology/cmorph.pyx
index c96039eb..4d491c3b 100644
--- a/skimage/morphology/cmorph.pyx
+++ b/skimage/morphology/cmorph.pyx
@@ -8,9 +8,9 @@ cimport numpy as np
 from libc.stdlib cimport malloc, free
 
 
-def _dilate(np.ndarray[np.uint8_t, ndim=2] image,
-            np.ndarray[np.uint8_t, ndim=2] selem,
-            np.ndarray[np.uint8_t, ndim=2] out=None,
+def _dilate(np.uint8_t[:, :] image,
+            np.uint8_t[:, :] selem,
+            np.uint8_t[:, :] out=None,
             char shift_x=0, char shift_y=0):
     """Return greyscale morphological dilation of an image.
 
@@ -52,12 +52,9 @@ def _dilate(np.ndarray[np.uint8_t, ndim=2] image,
     else:
         out = np.ascontiguousarray(out)
 
-    cdef np.uint8_t* out_data = <np.uint8_t*>out.data
-    cdef np.uint8_t* image_data = <np.uint8_t*>image.data
-
     cdef Py_ssize_t r, c, rr, cc, s, value, local_max
 
-    cdef Py_ssize_t selem_num = np.sum(selem != 0)
+    cdef Py_ssize_t selem_num = np.sum(np.asarray(selem) != 0)
     cdef Py_ssize_t* sr = <Py_ssize_t*>malloc(selem_num * sizeof(Py_ssize_t))
     cdef Py_ssize_t* sc = <Py_ssize_t*>malloc(selem_num * sizeof(Py_ssize_t))
 
@@ -76,21 +73,21 @@ def _dilate(np.ndarray[np.uint8_t, ndim=2] image,
                 rr = r + sr[s]
                 cc = c + sc[s]
                 if 0 <= rr < rows and 0 <= cc < cols:
-                    value = image_data[rr * cols + cc]
+                    value = image[rr, cc]
                     if value > local_max:
                         local_max = value
 
-            out_data[r * cols + c] = local_max
+            out[r, c] = local_max
 
     free(sr)
     free(sc)
 
-    return out
+    return np.asarray(out)
 
 
-def _erode(np.ndarray[np.uint8_t, ndim=2] image,
-           np.ndarray[np.uint8_t, ndim=2] selem,
-           np.ndarray[np.uint8_t, ndim=2] out=None,
+def _erode(np.uint8_t[:, :] image,
+           np.uint8_t[:, :] selem,
+           np.uint8_t[:, :] out=None,
            char shift_x=0, char shift_y=0):
     """Return greyscale morphological erosion of an image.
 
@@ -131,12 +128,9 @@ def _erode(np.ndarray[np.uint8_t, ndim=2] image,
     else:
         out = np.ascontiguousarray(out)
 
-    cdef np.uint8_t* out_data = <np.uint8_t*>out.data
-    cdef np.uint8_t* image_data = <np.uint8_t*>image.data
-
     cdef int r, c, rr, cc, s, value, local_min
 
-    cdef Py_ssize_t selem_num = np.sum(selem != 0)
+    cdef Py_ssize_t selem_num = np.sum(np.asarray(selem) != 0)
     cdef Py_ssize_t* sr = <Py_ssize_t*>malloc(selem_num * sizeof(Py_ssize_t))
     cdef Py_ssize_t* sc = <Py_ssize_t*>malloc(selem_num * sizeof(Py_ssize_t))
 
@@ -155,13 +149,13 @@ def _erode(np.ndarray[np.uint8_t, ndim=2] image,
                 rr = r + sr[s]
                 cc = c + sc[s]
                 if 0 <= rr < rows and 0 <= cc < cols:
-                    value = image_data[rr * cols + cc]
+                    value = image[rr, cc]
                     if value < local_min:
                         local_min = value
 
-            out_data[r * cols + c] = local_min
+            out[r, c] = local_min
 
     free(sr)
     free(sc)
 
-    return out
+    return np.asarray(out)
diff --git a/skimage/morphology/greyreconstruct.py b/skimage/morphology/greyreconstruct.py
index 049488f8..17455b55 100644
--- a/skimage/morphology/greyreconstruct.py
+++ b/skimage/morphology/greyreconstruct.py
@@ -159,7 +159,7 @@ def reconstruction(seed, mask, method='dilation', selem=None, offset=None):
     # Create a list of strides across the array to get the neighbors within
     # a flattened array
     value_stride = np.array(images.strides[1:]) / images.dtype.itemsize
-    image_stride = images.strides[0] / images.dtype.itemsize
+    image_stride = images.strides[0] // images.dtype.itemsize
     selem_mgrid = np.mgrid[[slice(-o, d - o)
                             for d, o in zip(selem.shape, offset)]]
     selem_offsets = selem_mgrid[:, selem].transpose()
diff --git a/skimage/segmentation/_felzenszwalb_cy.pyx b/skimage/segmentation/_felzenszwalb_cy.pyx
index 3babb86e..8590e17d 100644
--- a/skimage/segmentation/_felzenszwalb_cy.pyx
+++ b/skimage/segmentation/_felzenszwalb_cy.pyx
@@ -12,7 +12,8 @@ from skimage.morphology.ccomp cimport find_root, join_trees
 from ..util import img_as_float
 
 
-def _felzenszwalb_grey(image, double scale=1, sigma=0.8, Py_ssize_t min_size=20):
+def _felzenszwalb_grey(image, double scale=1, sigma=0.8,
+                       Py_ssize_t min_size=20):
     """Felzenszwalb's efficient graph based segmentation for a single channel.
 
     Produces an oversegmentation of a 2d image using a fast, minimum spanning
diff --git a/skimage/segmentation/_slic.pyx b/skimage/segmentation/_slic.pyx
index bfa6cd1a..5df747ad 100644
--- a/skimage/segmentation/_slic.pyx
+++ b/skimage/segmentation/_slic.pyx
@@ -2,15 +2,11 @@
 #cython: boundscheck=False
 #cython: nonecheck=False
 #cython: wraparound=False
-import collections as coll
 import numpy as np
-from time import time
 from scipy import ndimage
 
-cimport numpy as cnp
-
-from ..util import img_as_float, regular_grid
-from ..color import rgb2lab, gray2rgb
+from skimage.util import img_as_float, regular_grid
+from skimage.color import rgb2lab, gray2rgb
 
 
 def _slic_cython(double[:, :, :, ::1] image_zyx,
@@ -19,18 +15,18 @@ def _slic_cython(double[:, :, :, ::1] image_zyx,
                  double[:, ::1] means,
                  Py_ssize_t max_iter, Py_ssize_t n_segments):
     """Helper function for SLIC segmentation.
-    
+
     Parameters
     ----------
-    image_zyx : 4D np.ndarray of double, shape (Z, Y, X, 6)
+    image_zyx : 4D array of double, shape (Z, Y, X, 6)
         The image with embedded coordinates, that is, `image_zyx[i, j, k]` is
         `array([i, j, k, r, g, b])` or `array([i, j, k, L, a, b])`, depending
         on the colorspace.
-    nearest_mean : 3D np.ndarray of int, shape (Z, Y, X)
+    nearest_mean : 3D array of int, shape (Z, Y, X)
         The (initially empty) label field.
-    distance : 3D np.ndarray of double, shape (Z, Y, X)
+    distance : 3D array of double, shape (Z, Y, X)
         The (initially infinity) array of distances to the nearest centroid.
-    means : 2D np.ndarray of double, shape (n_segments, 6)
+    means : 2D array of double, shape (n_segments, 6)
         The centroids obtained by SLIC.
     max_iter : int
         The maximum number of k-means iterations.
@@ -39,7 +35,7 @@ def _slic_cython(double[:, :, :, ::1] image_zyx,
 
     Returns
     -------
-    nearest_mean : 3D np.ndarray of int, shape (Z, Y, X)
+    nearest_mean : 3D array of int, shape (Z, Y, X)
         The label field/superpixels found by SLIC.
     """
 
diff --git a/skimage/transform/_warps_cy.pyx b/skimage/transform/_warps_cy.pyx
index f51eb8d3..5dc2ffd0 100644
--- a/skimage/transform/_warps_cy.pyx
+++ b/skimage/transform/_warps_cy.pyx
@@ -83,10 +83,8 @@ def _warp_fast(cnp.ndarray image, cnp.ndarray H, output_shape=None,
 
     """
 
-    cdef cnp.ndarray[dtype=cnp.double_t, ndim=2, mode="c"] img = \
-         np.ascontiguousarray(image, dtype=np.double)
-    cdef cnp.ndarray[dtype=cnp.double_t, ndim=2, mode="c"] M = \
-         np.ascontiguousarray(H)
+    cdef double[:, ::1] img = np.ascontiguousarray(image, dtype=np.double)
+    cdef double[:, ::1] M = np.ascontiguousarray(H)
 
     if mode not in ('constant', 'wrap', 'reflect', 'nearest'):
         raise ValueError("Invalid mode specified.  Please use "
@@ -101,8 +99,7 @@ def _warp_fast(cnp.ndarray image, cnp.ndarray H, output_shape=None,
         out_r = output_shape[0]
         out_c = output_shape[1]
 
-    cdef cnp.ndarray[dtype=cnp.double_t, ndim=2] out = \
-         np.zeros((out_r, out_c), dtype=np.double)
+    cdef double[:, ::1] out = np.zeros((out_r, out_c), dtype=np.double)
 
     cdef Py_ssize_t tfr, tfc
     cdef double r, c
@@ -122,8 +119,8 @@ def _warp_fast(cnp.ndarray image, cnp.ndarray H, output_shape=None,
 
     for tfr in range(out_r):
         for tfc in range(out_c):
-            _matrix_transform(tfc, tfr, <double*>M.data, &c, &r)
-            out[tfr, tfc] = interp_func(<double*>img.data, rows, cols, r, c,
+            _matrix_transform(tfc, tfr, &M[0, 0], &c, &r)
+            out[tfr, tfc] = interp_func(&img[0, 0], rows, cols, r, c,
                                         mode_c, cval)
 
-    return out
+    return np.asarray(out)