3D Unet小簡介 — 實作(3)
終於來到到本系列的最後一篇了。回顧前幾篇，從理論到實作，其實有些地方還是講得不是很清楚。有些是礙於篇幅，有些則可能自己也還想得不夠仔細。不過，這些細節也可以之後再繼續寫系列文分享和討論。總之，本篇還是將照著預定的進度收尾。
本篇要討論的主要是兩個項目：model和train。
神經網路架構
有兩個檔案：unet和isensee2017。分別對應的是一般的3D Unet架構和改良過的模型。關於這兩個模型的差異請看這一篇。本篇不再花費篇幅描述架構差異，僅就程式碼的部分說明。
create_convolution_block
unet.py中最關鍵的是create_convolution_block。除了利用keras直接建構3D的CNN之外，還提到兩個normalization：batch normalization和instance normalization。因為函式有提到，所以稍作說明。想更深入了解可以搜尋台大電機李宏毅老師的投影片[1]。
if batch_normalization:
layer = BatchNormalization(axis=1)(layer)
elif instance_normalization:
try:
from keras_contrib.layers.normalization import
InstanceNormalizationnormalization中文翻做正規化，也就是「標準一致」的意思。什麼東西標準一致?特徵，或說梯度。未正規之前特徵的數值差異在不同軸向並不一致，可能在x軸是差1或10，在y軸是100。因此對應的權重和學習速率都有可能不同。如果可以透過正規化將差異調整成一致，學習會更有效率。
feature scaling是normalization的一種方法，其將每一筆資料減掉平均數再除以標準差。這個做法讓我想到中央極限定理(central limit theory)[2]，樣本數很大時其分佈會接近常態分佈(normal distribution)，也就是常態化的作法。啊…結果跟統計真的很有關聯。feature scaling用在神經網路，特別是深層的神經網路很有效，因為它可以將參數更新的幅度一致化，讓神經網路間的「溝通」變得更量化。但是，因為神經網路之間的參數總是在更動，要做到normalization是很不容易的。
batch normalization的精神就是以batch為單位，在每一層進入啟動函數前 — 對，強調之前 — 先做normalization。有些啟動函數例如sigmoid或是tanh，再曲線的兩端其實移動梯度很小，有點類似飽和的狀況。為了減少將輸出映射到這個區間的機率，可以用normalization。利用減掉平均數除以標準差的作法處理，然後再丟進啟動函數裡。平均值和標準差取決於本層網路的輸出值，所以和輸出值一樣，是會隨神經網路的反向傳播更動的。
經過normalization後的輸出(進入啟動函數前)的資料，其平均值會變成0，標準差則為1(根據中央極限定理)。如果不想把normalization的平均數和標準差都定為0/1的話，可以加一點參數。通常會讓輸出乘上一個值(γ)再加上一個值(β)，這是因為平均數受到β影響，而標準差受到γ影響。這兩個數是獨立的，類似神經網路的超參數(hyperparameter)，不會受到反向傳播的影響。
normalization的好處有：
1. 降低梯度消失的效應
2. 輔助某些啟動函數，強化其分辨能力
3. 降低權重等參數初始設定不佳帶來的負面影響
4. 有部分對抗過度配適(overfitting)的效果
那instance normalization又是什麼呢?batch normalization是以batch為單位來做校正，所以一個batch裡面必須有足夠的資料量，否則效果有限。但在某些情況下資料量是少的，或者資料量小反而訓練效果更好[3]，那就必須把normalization的操作方式限縮到以個別資料為單位來進行。所謂「個別資料」，指的是單通道(channel)、單樣本，以每個畫素為單位。
所以它和batch normalization的差異就是 — 想成一個批次只用一個個案來做batch normalization就行了。公式就不轉貼，可以到參考文章去看，基本上就差在「批次」這一項。此外，它也有β和γ。
照程式源碼看，keras似乎不支援instance normalization，而必須再另外下載擴充包安裝進去。不過tensorflow本身是有的，如果您的功力夠強，也可以直接用tensorflow呼叫。
最後，create_convolution_block啟動函數用的是ReLu。好像本身就叫不屬於易飽和的性質XD
get_up_convolution
get_up_convolution則是另一個關鍵。deconvolution和upsampling有什麼差異呢?這兩個做法其實最早都可以回溯到heat map的製作。當初為了要驗證卷積網路如何辨識圖片[4]，其根據為何?於是想到了「反置卷積網路」的做法，利用逆推讓便是圖片中的「熱區」現形。
deconvolution翻為「反卷積」，其實這個名詞不太正確。在操作上只是利用了矩陣的「轉置矩陣」(trasnpose)來運算而已，因此用transpose convolution這個說法更精確。所以它用特徵圖的一部分反推原本可能的分布是如何，而每一部分該反推多少就牽涉到其核(kernal)，一如我們在CNN中遇到的。因此這些圖是被「製造」出來，而非原先的圖形。
upsampling也是一樣，只是它用更單純的作法：用特徵圖中的值直接填滿該反推的區域。所以經由它生產出來的heat map看起來顆粒比較明顯，沒有滑順的邊界；但兩者在使用上是各有好壞。upsampling運算速度快；deconvolution則較精緻。
源碼提供兩者，可以自由切換。
unet_model_3d
來看主架構吧。順著註記看其實很清楚。
for layer_depth in range(depth):
layer1 = create_convolution_block(input_layer=current_layer,
n_filters=n_base_filters*(2**layer_depth),
batch_normalization=batch_normalization)
layer2 = create_convolution_block(input_layer=layer1,
n_filters=n_base_filters*(2**layer_depth)*2,
batch_normalization=batch_normalization)
if layer_depth < depth - 1:
current_layer = MaxPooling3D(pool_size=pool_size)
(layer2)
levels.append([layer1, layer2, current_layer])
else:
current_layer = layer2
levels.append([layer1, layer2])這一段是建立Unet的contracting path。每一階都會經過兩個卷積，然後只要還沒到最低階，就用池化層往下降。除此之外，把每一層也儲存起來，等一下concatenation要用到。
for layer_depth in range(depth-2, -1, -1):
up_convolution = get_up_convolution(pool_size=pool_size,
deconvolution=deconvolution,
n_filters=current_layer._keras_shape[1])(current_layer)
concat = concatenate([up_convolution, levels[layer_depth]
[1]], axis=1)
current_layer = create_convolution_block(n_filters=
levels[layer_depth][1]._keras_shape[1],
input_layer=concat,
batch_normalization=batch_normalization)
current_layer = create_convolution_block(n_filters=
levels[layer_depth][1]._keras_shape[1],
input_layer=current_layer,
batch_normalization=batch_normalization)這一段建立expanding path。從最底階開始，一層一層用deconvolution/upsampling升階，並和之前contracting path的對應同階層concatenation。直到回到最上階。
最後再用一個3D卷積輸出。損失用dice coefficient評估，這在前面已經講過了，這裡不再贅述。但要講一個optimzer:Adam。(81行)
model.compile(optimizer=Adam(lr=initial_learning_rate), .....所以「亞當」是什麼?有其他選擇嗎?例如「夏娃」….
其實這個「優化器」是keras裡面的其中一種設定，方便根據優化器內容調用各種與訓練相關的超參數，例如學習速率、動量、衰變等等。keras有內建數個優化器的模型，使用時可以先建立實例，設定好參數之後再於compile呼叫；也可以直接於compile呼叫指定優化器，但使用的就是keras裡預設的數值。順帶一提，裡面有各種模型名稱，但沒有夏娃XD想知道細節可以參閱keras的文件，這邊附上簡中版(找不到繁中)[5]。
跑完這個主程式，應該就可以建立一個3D Unet model實例了。
isensee2017的差異處
主要有幾個不同：
1. 模型架構
2. 預設取消instance normalization
3. 啟動函數改用leaky ReLu/sigmoid
模型的架構主要差在context model，以及只用upsampling、捨棄deconvolution。
若看其源碼，context module主要也是由兩層CNN組成，但多了一個dropout。表示視情況有些階可能會少一層convolution。以dropout的基本精神來看，應該是為了防止overfitting[6]。之後，在context module和前導的CNN層用summation相加，比較符合論文提出的原型。當然，和unet原型一樣，也要把每階的層儲存起來，等後面expanding path用。
expanding path的主軸是upsampling，但是：
1. 和前面contracting path的skip connection是用concatenation，而非summation。
2. 另外用segmentation layers把concatenation之後的層先經過一層CNN再儲存起來，準備後面圖形最右邊疊加時使用。
output_layer = None
for level_number in reversed(range(n_segmentation_levels)):
segmentation_layer = segmentation_layers[level_number]
if output_layer is None:
output_layer = segmentation_layer
else:
output_layer = Add()([output_layer, segmentation_layer]) if level_number > 0:
output_layer = UpSampling3D(size=(2, 2, 2)(output_layer)
這是最後一段，也就是上圖最右側的相加部分。也是用summation的方式完成。當然，因為兩階的維度不同，必須先upsampling之後將維度調整成一致，才能相加。
以上是兩個模型架構的程式碼介紹。接著來看看執行訓練用的檔案。
執行訓練吧！
根據原作者的建議，可以先使用isensee2017的模型訓練，之後再做調整。當然，也可以用unet3d的模型訓練，只要在import的地方調整就可以了。下面將以isensee2017的模型執行訓練。
利用config做成一組辭典，將必要的參數及所需的相關路徑指定好。執行的步驟從讀取資料開始。將資料按照之前所說的寫進hdf5檔，然後再讀取進data generator，最後進到神經網路訓練。訓練完成之後，本次訓練所拆分的training/validation data清單、訓練好的神經網路模型、以及資料hdf5檔都會被儲存致指定的路徑，以供後續訓練或測試使用。
在我的設定裡，是把patch shape關掉，labels設(1,)因為只有一種label(不知道該設(0,1)還是(1,))，mordality只有一組，所以nb_channels也是1。deconvolution設為true，batch size設很小(1，我知道對batch normalization來說這樣不理想，但應該不至於影響訓練?)。epoch只設5，全sample用第一批極少的10個資料先測試訓練能否順利完成，訓練成果則先不論。其餘跟patch相關的、data augmentation相關的，一律設為false。其他超參數按照預設值，並未調整。
訓練出了什麼問題?
遇到的問題是訓練跑完第一個epoch之後就卡在第二個，不再前進。同時跳出一個user warning：An input could not be retrieved. It could be because a worker has died.We do not have any information on the lost sample.
雖然不是error，可是訓練也沒辦法再進行下去了。只好把這段敘述拿去google。結果有各式各樣的答案：有人說是keras/tensorflow環境的問題，有人說是mask的形狀，有人說是記憶體，有人說是data generator的問題…不一而足。我試過的方法有：更換keras/tensorflow版本、重寫資料儲存檔案，從hdf5改成json、試著調用GPU，不過都沒效。現在不禁懷疑是data本身的問題，但至少array的維度和input都是吻合的。我不確定還有沒有其他檢查資料的方法，因為至少在進神經網路前的array我都有打開來看過，沒什麼異常。
從warning的敘述看起來，應該是某個階段的資料傳遞出現問題。究竟是第一輪的validation出現問題，還是第二輪的training出現問題?不得而知。因為在console中確實呈現了epoch 2/5這一行。但就沒有接下來的資訊了，有的只是上面那一行warning。試著從keras的source code尋找問題，但只能懷疑是在第一輪進到第二輪的地方出現錯誤，卻沒有進一步的頭緒。
還有沒有什麼其他可能呢?希望能有高手告訴我。
3D Unet的介紹也到這邊全部結束了。可惜最後沒有一個happy ending，但希望有人能幫我把結尾完成。如果有什麼想法或靈感，歡迎和我聯絡或提出。謝謝看到這邊的每位朋友。
Reference
[1] http://violin-tao.blogspot.com/2018/02/ml-batch-normalization.html
[2]http://www.math.nsysu.edu.tw/StatDemo/CentralLimitTheorem/CentralLimit.html
[3] https://zhuanlan.zhihu.com/p/56542480
[4] M.D.Zeiler, R.Fergus. Visualizing and Understanding Convolutional Networks, 2013.
[5] https://keras.io/zh/optimizers/
[6] A.Krizhevsky, I.Sutskever, G.E.Hinton. ImageNet Classification with Deep Convolutional.
