4. 影像的 NumPy 速成課程#

scikit-image 中的影像以 NumPy ndarray 表示。因此,許多常見的操作可以使用標準 NumPy 方法來操作陣列

>>> import skimage as ski
>>> camera = ski.data.camera()
>>> type(camera)
<type 'numpy.ndarray'>

注意

標記的類陣列資料類型,例如 pandas.DataFramexarray.DataArray,在 scikit-image 中不支援原生支援。但是,儲存在這些類型中的資料可以轉換為 numpy.ndarray 並帶有一些假設(請參閱 pandas.DataFrame.to_numpy()xarray.DataArray.data)。特別是,這些轉換會忽略取樣座標 (DataFrame.indexDataFrame.columnsDataArray.coords),這可能會導致資料表示錯誤,例如,當原始資料點不規則間隔時。

擷取影像的幾何形狀和像素數量

>>> camera.shape
(512, 512)
>>> camera.size
262144

擷取有關影像強度值的統計資訊

>>> camera.min(), camera.max()
(0, 255)
>>> camera.mean()
118.31400299072266

表示影像的 NumPy 陣列可以是不同的整數或浮點數值類型。請參閱 影像資料類型及其意義 以取得有關這些類型以及 scikit-image 如何處理它們的更多資訊。

4.1. NumPy 索引#

NumPy 索引可用於查看像素值和修改像素值

>>> # Get the value of the pixel at the 10th row and 20th column
>>> camera[10, 20]
153
>>> # Set to black the pixel at the 3rd row and 10th column
>>> camera[3, 10] = 0

請小心!在 NumPy 索引中,第一個維度 (camera.shape[0]) 對應於列,而第二個維度 (camera.shape[1]) 對應於行,原點 (camera[0, 0]) 位於左上角。這與矩陣/線性代數符號一致,但與笛卡爾 (x, y) 座標相反。 請參閱下方座標慣例以取得更多詳細資訊。

除了單個像素之外,還可以使用 NumPy 的不同索引功能來存取/修改整組像素的值。

切片

>>> # Set the first ten lines to "black" (0)
>>> camera[:10] = 0

遮罩 (使用布林遮罩索引)

>>> mask = camera < 87
>>> # Set to "white" (255) the pixels where mask is True
>>> camera[mask] = 255

花式索引 (使用索引集索引)

>>> import numpy as np
>>> inds_r = np.arange(len(camera))
>>> inds_c = 4 * inds_r % len(camera)
>>> camera[inds_r, inds_c] = 0

當您需要選擇一組像素來執行操作時,遮罩非常有用。遮罩可以是與影像形狀相同的任何布林陣列(或可廣播到影像形狀的形狀)。這可以用於定義感興趣的區域,例如磁碟

>>> nrows, ncols = camera.shape
>>> row, col = np.ogrid[:nrows, :ncols]
>>> cnt_row, cnt_col = nrows / 2, ncols / 2
>>> outer_disk_mask = ((row - cnt_row)**2 + (col - cnt_col)**2 >
...                    (nrows / 2)**2)
>>> camera[outer_disk_mask] = 0
../_images/sphx_glr_plot_camera_numpy_001.png

NumPy 中的布林運算可以用於定義更複雜的遮罩

>>> lower_half = row > cnt_row
>>> lower_half_disk = np.logical_and(lower_half, outer_disk_mask)
>>> camera = data.camera()
>>> camera[lower_half_disk] = 0

4.2. 彩色影像#

以上所有內容對彩色影像仍然成立。彩色影像是帶有額外尾部維度表示通道的 NumPy 陣列

>>> cat = ski.data.chelsea()
>>> type(cat)
<type 'numpy.ndarray'>
>>> cat.shape
(300, 451, 3)

這表示 cat 是一個 300 x 451 像素的影像,具有三個通道(紅色、綠色和藍色)。和以前一樣,我們可以取得和設定像素值

>>> cat[10, 20]
array([151, 129, 115], dtype=uint8)
>>> # Set the pixel at (50th row, 60th column) to "black"
>>> cat[50, 60] = 0
>>> # set the pixel at (50th row, 61st column) to "green"
>>> cat[50, 61] = [0, 255, 0]  # [red, green, blue]

我們也可以像上面處理灰階影像一樣,對 2D 多通道影像使用 2D 布林遮罩

(原始碼, png, hires.png, pdf)

../_images/numpy_images-1.png

在 2D 彩色影像上使用 2D 遮罩#

skimage.data 中包含的範例彩色影像會將通道儲存在最後一個軸上,儘管其他軟體可能會遵循不同的慣例。支援彩色影像的 scikit-image 函式有一個 channel_axis 引數,可用於指定陣列的哪個軸對應於通道。

4.3. 座標慣例#

因為 scikit-image 使用 NumPy 陣列來表示影像,因此座標慣例必須匹配。二維 (2D) 灰階影像(例如上面的 camera)由列和行索引(縮寫為 (row, col)(r, c)),左上角的最小元素為 (0, 0)。在程式庫的各個部分,您也會看到 rrcc 指的是列和行座標的清單。 我們將此慣例與 (x, y) 區分開來,(x, y) 通常表示標準笛卡爾座標,其中 x 是水平座標,y 是垂直座標,原點位於左下角(例如,Matplotlib 軸使用此慣例)。

對於多通道影像,任何維度(陣列軸)都可以用於顏色通道,並以 channelch 表示。在 scikit-image 0.19 之前,此通道維度始終是最後一個,但在目前的版本中,通道維度可以由 channel_axis 引數指定。需要多通道資料的函式預設為 channel_axis=-1。否則,函式預設為 channel_axis=None,表示假設沒有軸對應於通道。

最後,對於體積 (3D) 影像,例如影片、磁振造影 (MRI) 掃描、共聚焦顯微鏡等,我們將前導維度稱為 plane,縮寫為 plnp

這些慣例總結如下

scikit-image 中的維度名稱和順序慣例#

影像類型

座標

2D 灰階

(列, 行)

2D 多通道 (例如 RGB)

(列, 行, 通道)

3D 灰階

(平面, 列, 行)

3D 多通道

(平面, 列, 行, 通道)

請注意,ch 的位置由 channel_axis 引數控制。


scikit-image 中的許多函式可以直接對 3D 影像進行操作

>>> import numpy as np
>>> import scipy as sp
>>> import skimage as ski
>>> rng = np.random.default_rng()
>>> im3d = rng.random((100, 1000, 1000))
>>> seeds = sp.ndimage.label(im3d < 0.1)[0]
>>> ws = ski.segmentation.watershed(im3d, seeds)

但在許多情況下,第三個空間維度的解析度低於其他兩個維度。某些 scikit-image 函式提供 spacing 關鍵字引數以協助處理此類資料

>>> slics = ski.segmentation.slic(im3d, spacing=[5, 1, 1], channel_axis=None)

其他時候,必須以平面方式完成處理。當平面沿著前導維度堆疊時(符合我們的慣例),可以使用以下語法

>>> edges = np.empty_like(im3d)
>>> for pln, image in enumerate(im3d):
...     # Iterate over the leading dimension
...     edges[pln] = ski.filters.sobel(image)

4.4. 關於陣列維度順序的注意事項#

雖然軸的標籤看似任意,但它可能對運算速度產生顯著影響。這是因為現代處理器從記憶體中擷取資料時,從不會只擷取一個項目,而是擷取一整塊相鄰的項目(一種稱為預取的操作)。因此,在記憶體中彼此相鄰的元素的處理速度,會比分散的元素處理速度快,即使運算次數相同也是如此。

>>> def in_order_multiply(arr, scalar):
...     for plane in list(range(arr.shape[0])):
...         arr[plane, :, :] *= scalar
...
>>> def out_of_order_multiply(arr, scalar):
...     for plane in list(range(arr.shape[2])):
...         arr[:, :, plane] *= scalar
...
>>> import time
>>> rng = np.random.default_rng()
>>> im3d = rng.random((100, 1024, 1024))
>>> t0 = time.time(); x = in_order_multiply(im3d, 5); t1 = time.time()
>>> print("%.2f seconds" % (t1 - t0))  
0.14 seconds
>>> s0 = time.time(); x = out_of_order_multiply(im3d, 5); s1 = time.time()
>>> print("%.2f seconds" % (s1 - s0))  
1.18 seconds
>>> print("Speedup: %.1fx" % ((s1 - s0) / (t1 - t0)))  
Speedup: 8.6x

當最末端/最右側維度變得更大時,加速效果會更加顯著。在開發演算法時,值得思考資料局部性。特別是,scikit-image 預設使用 C 連續陣列。當使用巢狀迴圈時,陣列的最末端/最右側維度應該位於計算的最內層迴圈中。在上面的例子中,*= numpy 運算符會迭代所有剩餘的維度。

4.5. 關於時間維度的說明#

雖然 scikit-image 目前沒有提供專門處理時變 3D 資料的函式,但它與 NumPy 陣列的相容性使我們可以很自然地處理形狀為 (t, pln, row, col, ch) 的 5D 陣列。

>>> for timepoint in image5d:  
...     # Each timepoint is a 3D multichannel image
...     do_something_with(timepoint)

然後我們可以將上表補充如下

scikit-image 中維度名稱和順序的增補#

影像類型

坐標

2D 彩色影片

(t, row, col, ch)

3D 彩色影片

(t, pln, row, col, ch)