Sorry, but I still don't understand why peak_local_max would not work with N-D arrays... I just tested:

```
In [1]: from skimage.feature import peak_local_max

In [2]: import numpy as np

In [3]: a = np.zeros((20, 20, 20))

In [4]: a[10, 10, 10] = 1

In [5]: peak_local_max(a, min_distance=1)
Out[5]: array([[10, 10, 10]])
```

Though, when looking at the code, there might be a bug in the logic of `exclude_border`, shouldn't we set the values to -INF in lines 136, 137?