注意
前往結尾下載完整範例程式碼。或透過 Binder 在您的瀏覽器中執行此範例
估計 3D 顯微影像中的異向性#
在本教學中，我們計算 3D 影像的結構張量。有關 3D 影像處理的一般介紹，請參閱 探索 3D 影像 (細胞)。我們在此使用的資料是從共軛焦螢光顯微鏡獲得的腎臟組織影像中採樣的 (更多詳細資訊請參閱 [1] 在 kidney-tissue-fluorescence.tif 下)。
import matplotlib.pyplot as plt
import numpy as np
import plotly.express as px
import plotly.io
from skimage import data, feature
載入影像#
此生物醫學影像可透過 scikit-image 的資料登錄取得。
data = data.kidney()
我們的 3D 多通道影像的形狀和大小究竟是什麼？
print(f'number of dimensions: {data.ndim}')
print(f'shape: {data.shape}')
print(f'dtype: {data.dtype}')
number of dimensions: 4
shape: (16, 512, 512, 3)
dtype: uint16
在本教學中，我們僅考慮第二個色彩通道，這為我們留下了一個 3D 單通道影像。值的範圍是多少？
range: (68, 4095)
讓我們視覺化我們的 3D 影像的中間切片。
讓我們選擇一個顯示相對 X-Y 等向性的特定區域。相反地，沿著 Z 的梯度非常不同 (而且，就此而言，很弱)。
sample = data[5:13, 380:410, 370:400, 1]
step = 3
cols = sample.shape[0] // step + 1
_, axes = plt.subplots(nrows=1, ncols=cols, figsize=(16, 8))
for it, (ax, image) in enumerate(zip(axes.flatten(), sample[::step])):
ax.imshow(image, cmap='gray', vmin=v_min, vmax=v_max)
ax.set_title(f'Plane = {5 + it * step}')
ax.set_xticks([])
ax.set_yticks([])
若要以 3D 檢視樣本資料，請執行下列程式碼
import plotly.graph_objects as go
(n_Z, n_Y, n_X) = sample.shape
Z, Y, X = np.mgrid[:n_Z, :n_Y, :n_X]
fig = go.Figure(
data=go.Volume(
x=X.flatten(),
y=Y.flatten(),
z=Z.flatten(),
value=sample.flatten(),
opacity=0.5,
slices_z=dict(show=True, locations=[4])
)
)
fig.show()
計算結構張量#
讓我們視覺化樣本資料的底部切片，並確定強烈變化的典型大小。我們應將此大小用作視窗函數的「寬度」。
關於最亮的區域 (即，在列 ~ 22 和欄 ~ 17 處)，我們可以看到跨欄 (分別為跨列) 的 2 或 3 (分別為 1 或 2) 個像素的變化 (因此，有強烈的梯度)。因此，我們可能選擇，例如，sigma = 1.5 作為視窗函數。或者，我們可以按每個軸傳遞 sigma，例如，sigma = (1, 2, 3)。請注意，沿著第一個 (Z，平面) 軸的大小 1 聽起來很合理，因為後者的大小為 8 (13 - 5)。在 X-Z 或 Y-Z 平面中檢視切片確認這是合理的。
然後我們可以計算結構張量的特徵值。
(3, 8, 30, 30)
最大特徵值在哪裡？
coords = np.unravel_index(eigen.argmax(), eigen.shape)
assert coords[0] == 0 # by definition
coords
(np.int64(0), np.int64(1), np.int64(22), np.int64(16))
注意
讀者可以檢查此結果 (座標 (plane, row, column) = coords[1:]) 對於不同的 sigma 有多穩健。
讓我們檢視在發現最大特徵值 (即，Z = coords[1]) 的 X-Y 平面中特徵值的空間分佈。