本筆記本通過一個示例說明如何在非rgb數據上使用DenseCRFs。同時，它將解釋基本概念並通過一個示例進行演示，因此即使您正在處理RGB數據，它也可能是有用的，不過也請查看PyDenseCRF's README !
在jupyter notebook上運行https://github.com/lucasb-eyer/pydensecrf/tree/master/examples/Non RGB Example.ipynb
Basic setup
先導入：
import pydensecrf.densecrf as dcrf from pydensecrf.utils import unary_from_softmax, create_pairwise_bilateral
設置畫圖風格：
import numpy as np import matplotlib.pyplot as plt %matplotlib inline plt.rcParams['image.interpolation'] = 'nearest' #設置插入風格 plt.rcParams['image.cmap'] = 'gray' #設置顏色風格
Unary Potential一元勢
一元勢由每個像素的類概率組成。這可以來自任何類型的模型，如隨機森林或深度神經網絡的softmax。
1）Create unary potential
補充：scipy.stats.multivariate_normal的使用
1.得到x,y的取值
from scipy.stats import multivariate_normal H, W, NLABELS = 400, 512, 2 # This creates a gaussian blob... a = np.mgrid[0:H, 0:W] #(2,400,512) #得到x,y的取值 a,a.shape
返回：
(array([[[ 0, 0, 0, ..., 0, 0, 0], [ 1, 1, 1, ..., 1, 1, 1], [ 2, 2, 2, ..., 2, 2, 2], ..., [397, 397, 397, ..., 397, 397, 397], [398, 398, 398, ..., 398, 398, 398], [399, 399, 399, ..., 399, 399, 399]], [[ 0, 1, 2, ..., 509, 510, 511], [ 0, 1, 2, ..., 509, 510, 511], [ 0, 1, 2, ..., 509, 510, 511], ..., [ 0, 1, 2, ..., 509, 510, 511], [ 0, 1, 2, ..., 509, 510, 511], [ 0, 1, 2, ..., 509, 510, 511]]]), (2, 400, 512))
2.得到點(x,y)的值
pos = np.stack(a, axis=2) #得到(x,y)點的值 pos,pos.shape
返回：
(array([[[ 0, 0], [ 0, 1], [ 0, 2], ..., [ 0, 509], [ 0, 510], [ 0, 511]], [[ 1, 0], [ 1, 1], [ 1, 2], ..., [ 1, 509], [ 1, 510], [ 1, 511]], [[ 2, 0], [ 2, 1], [ 2, 2], ..., [ 2, 509], [ 2, 510], [ 2, 511]], ..., [[397, 0], [397, 1], [397, 2], ..., [397, 509], [397, 510], [397, 511]], [[398, 0], [398, 1], [398, 2], ..., [398, 509], [398, 510], [398, 511]], [[399, 0], [399, 1], [399, 2], ..., [399, 509], [399, 510], [399, 511]]]), (400, 512, 2))
3.得到概率值：
#mean=[200, 256], cov=12800,12800的單位乘作為協方差矩陣 rv = multivariate_normal([H//2, W//2], (H//4)*(W//4)) probs = rv.pdf(pos) probs,probs.shape
返回：
(array([[2.01479637e-07, 2.05541767e-07, 2.09669414e-07, ..., 2.13863243e-07, 2.09669414e-07, 2.05541767e-07], [2.04644486e-07, 2.08770423e-07, 2.12962908e-07, ..., 2.17222613e-07, 2.12962908e-07, 2.08770423e-07], [2.07842810e-07, 2.12033230e-07, 2.16291237e-07, ..., 2.20617517e-07, 2.16291237e-07, 2.12033230e-07], ..., [2.11074628e-07, 2.15330207e-07, 2.19654423e-07, ..., 2.24047973e-07, 2.19654423e-07, 2.15330207e-07], [2.07842810e-07, 2.12033230e-07, 2.16291237e-07, ..., 2.20617517e-07, 2.16291237e-07, 2.12033230e-07], [2.04644486e-07, 2.08770423e-07, 2.12962908e-07, ..., 2.17222613e-07, 2.12962908e-07, 2.08770423e-07]]), (400, 512))
4.將概率值進行處理集中到[0.4,0.6]區間中
先將值集中到[0,1]區間中
# ...which we project into the range [0.4, 0.6] probs = (probs-probs.min()) / (probs.max()-probs.min()) probs, probs.shape
返回：
(array([[0. , 0.00033208, 0.00066951, ..., 0.00101235, 0.00066951, 0.00033208], [0.00025872, 0.00059602, 0.00093875, ..., 0.00128698, 0.00093875, 0.00059602], [0.00052019, 0.00086275, 0.00121084, ..., 0.00156451, 0.00121084, 0.00086275], ..., [0.00078439, 0.00113228, 0.00148578, ..., 0.00184495, 0.00148578, 0.00113228], [0.00052019, 0.00086275, 0.00121084, ..., 0.00156451, 0.00121084, 0.00086275], [0.00025872, 0.00059602, 0.00093875, ..., 0.00128698, 0.00093875, 0.00059602]]), (400, 512))
再將其集中到[0.4,0.6]區間中，0.5 + 0.2 * (probs-0.5)等價於0.4+0.2*probs，而此時的probs為[0,1]區間中的值：
probs = 0.5 + 0.2 * (probs-0.5) #此時probs.shape為(400, 512) probs, probs.shape #將最終的概率結果的值集中在[0.4, 0.6]區間中
返回：
(array([[0.4 , 0.40006642, 0.4001339 , ..., 0.40020247, 0.4001339 , 0.40006642], [0.40005174, 0.4001192 , 0.40018775, ..., 0.4002574 , 0.40018775, 0.4001192 ], [0.40010404, 0.40017255, 0.40024217, ..., 0.4003129 , 0.40024217, 0.40017255], ..., [0.40015688, 0.40022646, 0.40029716, ..., 0.40036899, 0.40029716, 0.40022646], [0.40010404, 0.40017255, 0.40024217, ..., 0.4003129 , 0.40024217, 0.40017255], [0.40005174, 0.4001192 , 0.40018775, ..., 0.4002574 , 0.40018775, 0.4001192 ]]), (400, 512))
5.一份數據已經准備好了，要復制生成一個效果相反的數據：
# 第一個維度需要等於類的數量，這里設置兩個類，為2 # 讓我們有一個"foreground"和一個"background"類 #因此，復制高斯blob，但倒置它去創建與“background”類的概率是相反的“foreground”類。 #np.tile(a,(m,n)):即是把a數組里面的元素復制n次放進一個數組c中，然后再把數組c復制m次放進一個數組b中 #np.tile(a,(2,3)) #將a數組重復3次形成2行的數組 #np.newaxis為多維數組增加一個軸 #probs[np.newaxis,:,:].shape為(1, 400, 512) probs = np.tile(probs[np.newaxis,:,:],(2,1,1))#此時probs.shape為(2, 400, 512) probs[1,:,:] = 1 - probs[0,:,:] #得到相反的值
probs,probs.shape
返回：
(array([[[0.4 , 0.40006642, 0.4001339 , ..., 0.40020247, 0.4001339 , 0.40006642], [0.40005174, 0.4001192 , 0.40018775, ..., 0.4002574 , 0.40018775, 0.4001192 ], [0.40010404, 0.40017255, 0.40024217, ..., 0.4003129 , 0.40024217, 0.40017255], ..., [0.40015688, 0.40022646, 0.40029716, ..., 0.40036899, 0.40029716, 0.40022646], [0.40010404, 0.40017255, 0.40024217, ..., 0.4003129 , 0.40024217, 0.40017255], [0.40005174, 0.4001192 , 0.40018775, ..., 0.4002574 , 0.40018775, 0.4001192 ]], [[0.6 , 0.59993358, 0.5998661 , ..., 0.59979753, 0.5998661 , 0.59993358], [0.59994826, 0.5998808 , 0.59981225, ..., 0.5997426 , 0.59981225, 0.5998808 ], [0.59989596, 0.59982745, 0.59975783, ..., 0.5996871 , 0.59975783, 0.59982745], ..., [0.59984312, 0.59977354, 0.59970284, ..., 0.59963101, 0.59970284, 0.59977354], [0.59989596, 0.59982745, 0.59975783, ..., 0.5996871 , 0.59975783, 0.59982745], [0.59994826, 0.5998808 , 0.59981225, ..., 0.5997426 , 0.59981225, 0.5998808 ]]]), (2, 400, 512))
6.最終的代碼為：
from scipy.stats import multivariate_normal H, W, NLABELS = 400, 512, 2 # This creates a gaussian blob... pos = np.stack(np.mgrid[0:H, 0:W], axis=2) rv = multivariate_normal([H//2, W//2], (H//4)*(W//4)) probs = rv.pdf(pos) # ...which we project into the range [0.4, 0.6] probs = (probs-probs.min()) / (probs.max()-probs.min()) probs = 0.5 + 0.2 * (probs-0.5) #此時probs.shape為(400, 512) # 第一個維度需要等於類的數量 # 讓我們有一個"foreground"和一個"background"類 #因此，復制高斯blob，但倒置它去創建與“background”類的概率是相反的“foreground”類。 #np.tile(a,(m,n)):即是把a數組里面的元素復制n次放進一個數組c中，然后再把數組c復制m次放進一個數組b中 #np.tile(a,(2,3)) #將a數組重復3次形成2行的數組 #np.newaxis為多維數組增加一個軸 #probs[np.newaxis,:,:].shape為(1, 400, 512) probs = np.tile(probs[np.newaxis,:,:],(2,1,1))#此時probs.shape為(2, 400, 512) probs[1,:,:] = 1 - probs[0,:,:] #得到相反的值 # Let's have a look: plt.figure(figsize=(15,5)) plt.subplot(1,2,1); plt.imshow(probs[0,:,:]); plt.title('Foreground probability'); plt.axis('off'); plt.colorbar(); plt.subplot(1,2,2); plt.imshow(probs[1,:,:]); plt.title('Background probability'); plt.axis('off'); plt.colorbar();
返回：
2）Run inference with unary potential 使用一元勢運行推理
我們已經可以運行一個只有一元勢的DenseCRF。這不是一個好主意，但我們可以做到:
1.unary_from_softmax函數的作用
將softmax類概率轉換為一元勢(每個節點的NLL)。
即我們之前先對圖片使用訓練好的網絡預測得到最終經過softmax函數得到的分類結果，這里需要將這個結果轉成一元勢
參數
- sm: numpy.array ,第一個維度是類的softmax的輸出，其他所有維度都是flattend。這意味着“sm.shape[0] == n_classes”。
- scale: float,softmax輸出的確定性(默認為None)，需要值在(0,1]。如果不為None，則softmax輸出被縮放到從[0,scale]概率的范圍。
- clip: float,將概率裁剪到的最小值。這是因為一元函數是概率的負對數，而log(0) = inf，所以我們需要把0概率裁剪成正的值。
在這里因為scale=None,clip=None,所以這個函數的作用其實只進行了下面的操作：
-np.log(sm).reshape([num_cls, -1]).astype(np.float32)
對值取負對數，並將結果變為(2, 204800)大小，進行了flatten操作，將圖片壓扁了：
# Inference without pair-wise terms #設置固定的一元勢 from pydensecrf.utils import unary_from_softmax U = unary_from_softmax(probs) # note: num classes is first dim，類別在第一維 U,U.shape
返回：
(array([[0.91629076, 0.9161247 , 0.915956 , ..., 0.91564745, 0.9158215 , 0.9159928 ], [0.51082563, 0.5109363 , 0.5110488 , ..., 0.5112547 , 0.5111386 , 0.5110243 ]], dtype=float32), (2, 204800))
此時的值就是一元勢了
d = dcrf.DenseCRF2D(W, H, NLABELS) d.setUnaryEnergy(U) # Run inference for 10 iterations，10次迭代 Q_unary = d.inference(10) # Q現在是近似后驗，我們可以用argmax得到一個MAP估計值。 #np.argmax得到Q_unary列方向上的最大值的索引 map_soln_unary = np.argmax(Q_unary, axis=0)#此時形式為(204800,)
map_soln_unary,map_soln_unary.shape
返回：
(array([1, 1, 1, ..., 1, 1, 1]), (204800,))
# 不幸的是，DenseCRF把所有東西都壓扁了，所以把它恢復到圖片形式。 map_soln_unary = map_soln_unary.reshape((H,W)) #重新變為(400, 512)) map_soln_unary,map_soln_unary.shape
返回：
(array([[1, 1, 1, ..., 1, 1, 1], [1, 1, 1, ..., 1, 1, 1], [1, 1, 1, ..., 1, 1, 1], ..., [1, 1, 1, ..., 1, 1, 1], [1, 1, 1, ..., 1, 1, 1], [1, 1, 1, ..., 1, 1, 1]]), (400, 512))
總的代碼為：
# Inference without pair-wise terms U = unary_from_softmax(probs) # note: num classes is first dim d = dcrf.DenseCRF2D(W, H, NLABELS) d.setUnaryEnergy(U) # Run inference for 10 iterations，10次迭代 Q_unary = d.inference(10) # Q現在是近似后驗，我們可以用argmax得到一個MAP估計值。 #np.argmax得到Q_unary列方向上的最大值的索引 map_soln_unary = np.argmax(Q_unary, axis=0)#此時形式為(204800,) # 不幸的是，DenseCRF把所有東西都壓扁了，所以把它恢復到圖片形式。 map_soln_unary = map_soln_unary.reshape((H,W))#重新變為(400, 512)) # And let's have a look. plt.imshow(map_soln_unary); plt.axis('off'); plt.title('MAP Solution without pairwise terms');
返回：
因為設置了顏色為gray，實際結果為：
Pairwise terms二元
DenseCRFs的全部意義在於使用某種形式的內容來平滑預測。這是通過“成對”術語來實現的，“成對”術語編碼元素之間的關系。
1）Add (non-RGB) pairwise term
例如，在圖像處理中，一個流行的成對關系是“雙邊”(bilateral)關系，它大致表示顏色相似或位置相似的像素可能屬於同一個類。
NCHAN=1 #通道為1，non-RGB #創造簡單的雙邊圖象。 #注意，我們將通道維度放在最后，但我們也可以將它作為第一個維度，只要將chdim參數進一步更改為0。 img = np.zeros((H,W,NCHAN), np.uint8) img[H//3:2*H//3,W//4:3*W//4,:] = 1 #使得[133:266,128:384,:]內的值為1，下面畫為白色 plt.imshow(img[:,:,0]); plt.title('Bilateral image'); plt.axis('off'); plt.colorbar();
返回：
2.create_pairwise_bilateral函數的作用
創造成對雙邊勢的Util函數。這適用於所有的圖像尺寸。對於2D情況，與“densecrf2 . addpairwisebilateral”函數效果相同。
參數：
- sdims: list or tuple，每個維度的比例因子。這在DenseCRF2D.addPairwiseBilateral中稱為“sxy”。
- schan: list or tuple，圖像中每個通道的比例因子。這在DenseCRF2D.addPairwiseBilateral中稱為“srgb”。
- img: numpy.array，輸入的圖像
- chdim: int, 可選。這指定了通道維度在圖像中的位置。例如，' chdim=2 '用於大小為(240,300,3)的RGB圖像，說明其通道3值放在維度2上(維度從0開始)。如果圖像沒有通道尺寸(例如只有一個通道)，則使用' chdim=-1 '。
#從上圖中創建成對的雙邊術語。 #這兩個' s{dims,chan} '參數是模型超參數，分別定義了位置的強度和圖像內容的雙邊參數。 #將通道為1的圖變為了通道為3的圖 pairwise_energy = create_pairwise_bilateral(sdims=(10,10), schan=(0.01,), img=img, chdim=2) #其shape為(3, 204800) #pairwise_energy現在包含的維度和DenseCRF的特性一樣多，在本例中是3:(x,y,channel1) img_en = pairwise_energy.reshape((-1, H, W)) # Reshape just for plotting,為了畫圖 #下面將三個通道對應的圖畫出來 plt.figure(figsize=(15,5)) plt.subplot(1,3,1); plt.imshow(img_en[0]); plt.title('Pairwise bilateral [x]'); plt.axis('off'); plt.colorbar(); plt.subplot(1,3,2); plt.imshow(img_en[1]); plt.title('Pairwise bilateral [y]'); plt.axis('off'); plt.colorbar(); plt.subplot(1,3,3); plt.imshow(img_en[2]); plt.title('Pairwise bilateral [c]'); plt.axis('off'); plt.colorbar();
返回：
2）Run inference of complete DenseCRF
現在我們可以創建一個包含一元勢和成對勢的dense CRF ，並對其進行推理以得到最終結果。
d = dcrf.DenseCRF2D(W, H, NLABELS) d.setUnaryEnergy(U) #設置一元勢 #添加二元勢 d.addPairwiseEnergy(pairwise_energy, compat=10) # `compat` is the "strength" of this potential. #這一次，讓我們自己分步驟進行推理，這樣我們就可以查看中間解，並監視KL-divergence，它表明我們收斂得有多好。 #PyDenseCRF還要求我們跟蹤計算所需的兩個臨時緩沖區，tmp1, tmp2。 Q, tmp1, tmp2 = d.startInference() for _ in range(5): #迭代5次后，查看結果 d.stepInference(Q, tmp1, tmp2) kl1 = d.klDivergence(Q) / (H*W) map_soln1 = np.argmax(Q, axis=0).reshape((H,W)) for _ in range(20): d.stepInference(Q, tmp1, tmp2) kl2 = d.klDivergence(Q) / (H*W) map_soln2 = np.argmax(Q, axis=0).reshape((H,W)) for _ in range(50): d.stepInference(Q, tmp1, tmp2) kl3 = d.klDivergence(Q) / (H*W) map_soln3 = np.argmax(Q, axis=0).reshape((H,W)) img_en = pairwise_energy.reshape((-1, H, W)) # Reshape just for plotting plt.figure(figsize=(15,5)) plt.subplot(1,3,1); plt.imshow(map_soln1); plt.title('MAP Solution with DenseCRF\n(5 steps, KL={:.2f})'.format(kl1)); plt.axis('off'); plt.subplot(1,3,2); plt.imshow(map_soln2); plt.title('MAP Solution with DenseCRF\n(20 steps, KL={:.2f})'.format(kl2)); plt.axis('off'); plt.subplot(1,3,3); plt.imshow(map_soln3); plt.title('MAP Solution with DenseCRF\n(50 steps, KL={:.2f})'.format(kl3)); plt.axis('off');
返回：
這個不太好看出效果，下面：
可見50次迭代后效果很好