---
# System prepended metadata

title: 明龍的專題生們 Meeting專用

---

---
title: 明龍的專題生們 Meeting專用
---
### [專題簡介鏈接](https://docs.google.com/document/d/1CJc4wW1l27xfuhmg2Jr0-E8wjnMQnpgIC6-kBB91xdQ/edit)
# Dental Segmentation in Cone-beam Computed Tomography Images
## Preprocessing

|                Input image                 |                Thresholding                |              Gamma correction              |
|:------------------------------------------:|:------------------------------------------:|:------------------------------------------:|
| ![](https://i.imgur.com/dSPRXYO.png =150x) | ![](https://i.imgur.com/22Cdz2l.png =150x) | ![](https://i.imgur.com/Huy4w4n.png =150x) |

- Thresholding: 去除背景的噪點以取得清晰的影像
- Gamma correction: 使得牙齒以及其他非牙齒部分亮度有較明顯對比


|                Whole image                 |MIP VALUE|               Image cropped                |
|:------------------------------------------:|:-:|:------------------------------------------:|
| ![](https://i.imgur.com/IpfZGcl.png =200x) |![](https://i.imgur.com/ZY3tdq5.png)|![](https://i.imgur.com/u7z57Tx.png =200x) |

- 實作方式: 先對影像作 Max Intensity Projection，取得影像的樣子(上左圖)，而後留下整張影像中亮度最亮的前10%，去除離群值後計算位置的平均作為圖片的中心只留下前後分別 150 slices

## Augmentation
- 因為我們擁有的 label data 過少，因此我們必須要利用旋轉或者鏡射的方式使得資料更為豐富

|       |                  Original                  |              Horizontal Flip               |                   Rotate                   | Zoom |
|:-----:|:------------------------------------------:|:------------------------------------------:|:------------------------------------------:|:----:|
| Image | ![](https://i.imgur.com/Huy4w4n.png =150x) | ![](https://i.imgur.com/curBedj.png =150x) | ![](https://i.imgur.com/9at9k7G.png =150x) |![](https://i.imgur.com/CKBqXPZ.png =150x)|
| Label | ![](https://i.imgur.com/I60s70w.png =150x) | ![](https://i.imgur.com/lnGbKgZ.png =150x) | ![](https://i.imgur.com/xgWvKP1.png =150x) |![](https://i.imgur.com/XsqzYUO.png =150x)|

## Model structure
![](https://i.imgur.com/00koiHk.png)

- 在 Unet 的 structure 上，我們加入 batch normalization 使得 hidden layer 比較穩定，並且將影像取 canny 的結果與 unet 結果連接之後作 convolution，以取得更好的結果
- loss function: 我們使用 dice_loss 以確保在兩個 class 數量相差懸殊還是可以有較好的訓練結果。
    $dice\space loss =1 - \dfrac{2\times(y_{true}\cap y_{pred})}{y_{true} + y_{pred}}$
    
## Post processing

* 我們將預測出來的300張slice疊成3d影像，並且再分別對x, y座標分割，各得到512張平面影像，總計兩份。
* 使用morphology的opening運算，依序照這兩個方向處理
* 然後再把圖片中的離群值去掉
* 但目前還有問題，其中一個就是kernel大小，還有不清楚morpholy的可行性。<br />從圖中可見，進行處理後雖然可以清除雜訊，卻也有一些不能刪掉的部分消失了。而有一些雜訊卻沒有消失。

|                before                |                after                 |
|:------------------------------------:|:------------------------------------:|
| ![](https://i.imgur.com/bKLcei8.png) | ![](https://i.imgur.com/ytG4WLN.png) |
|  |  |

## Results
STL檔
https://www.embodi3d.com/files/file/47177-final0059-stl-file-processed/
https://www.embodi3d.com/files/file/47154-final0005-stl-file-processed/
https://www.embodi3d.com/files/file/47162-final0011-stl-file-processed/

# 其他嘗試

## 1. 將canny edge 從一開始就和原圖一起訓練， 生成結果
![](https://i.imgur.com/1kRkgoK.png)
![](https://i.imgur.com/Nr0RYua.png)

![](https://i.imgur.com/cNsYoik.png)

## Epoch = 30
![](https://i.imgur.com/qPMAJwF.png)![](https://i.imgur.com/eSqSlKt.png)

## Epoch = 10
![](https://i.imgur.com/LoNs08N.png)
![](https://i.imgur.com/1ML4c4i.png)
### 結論： canny edge 放到前面會造成牙根分割效果降低，所以依然采用canny edge 放在後面的效果

## 2. HU window 延伸處理



| 原圖 | 轉HU | 二值化 |
| -------- | -------- | -------- |
|![](https://i.imgur.com/Kn8ouOa.png)|![](https://i.imgur.com/U8gr0BQ.png)| ![](https://i.imgur.com/wIvRhh3.png)|
|連通區域標記|去除碎片|腐蝕操作(Erosion)|
![](https://i.imgur.com/pjwvO3K.png)|![](https://i.imgur.com/qD88J1m.png)|![](https://i.imgur.com/A8KxAVB.png)|
|閉合運算|孔洞填充|對應原圖|
![](https://i.imgur.com/Tys6Klm.png)|![](https://i.imgur.com/wO4pOPL.png)|![](https://i.imgur.com/j5arKPi.png)|

二值化： 去HU值低於800的為0, 高於800的為1, 700以上為骨頭,800分割效果較好
連通區域標記： 把相連一起的點視爲同個區域
去除碎片：保留16個最大聯通區域(上下顎各16顆牙齒)，其餘的設爲0
腐蝕操作：將牙齒與牙槽骨相連處斷開
閉合運算：去掉牙槽骨
孔洞填充：填滿牙齒内小孔
對應原圖：投影到原圖上

|原前處理後結果|經上述前處理後結果|
| - | - |
|![](https://i.imgur.com/7klzqjt.png)|![](https://i.imgur.com/8yhotzW.png)|

* ![](https://i.imgur.com/kq9oH0o.png)
* ![](https://i.imgur.com/Mrij202.png)

* 參考網站説明： [常見醫療掃描影象處理步驟](https://www.itread01.com/content/1544094552.html)

### 結論： 認爲分割效果不錯，另作爲feature之一，同canny input一同加到訓練模型裏

### 延伸：如果合并之前的前處理（threshold & gmm) 呢？
* 因爲threshold 和 gmm 會更變原圖的值，造成計算HU時會有誤，目前沒有時間再去校正，所以否決這個辦法
* ![](https://i.imgur.com/7G2plzK.png)

### 延伸：如果只使用切割後的結果代替原來的input呢？
* 下面是切割結果
* ![](https://i.imgur.com/Tk69fGU.png)
* ![](https://i.imgur.com/mNbroQL.png)
* 推測：HU 值加上切割可能失去原圖的一些資訊，造成訓練效果較不佳


# transfer learning的後處理
case57 slide 177
case93 slide 141
|  Input image  |  預測結果  |    後處理結果     |
|--------------|------------|------------|
| ![](https://i.imgur.com/6TLiLws.png) | ![](https://i.imgur.com/EjO9vri.png) | ![](https://i.imgur.com/uKLiQ4h.png) |
| ![](https://i.imgur.com/C6dzBmx.png) |![](https://i.imgur.com/w8M1hcv.png) | ![](https://i.imgur.com/qU0ZTRZ.png) |
經過transfer learing後，只剩下少數case有明顯雜訊。所以後處理的重要性有些降低，主要拿來去除明顯離群值和一些雜訊，但也因此得以自動化

# 有無transfer的差別
case 57_177
case 57_166
case 57_161
可以看出有transfer後，需要後處理的部分減少很多，預測變得更精準
|  Input image  |  無transfer預測結果  |    無transfer後處理結果     |  transfer預測結果  |    transfer後處理結果     |
|--------------|------------|------------|------------|------------|
| ![](https://i.imgur.com/KKGSHZ1.png) | ![](https://i.imgur.com/VHUPCbl.png) | ![](https://i.imgur.com/u9AvevW.png) | ![](https://i.imgur.com/OO1LMj9.png) | ![](https://i.imgur.com/QKFpqae.png) |
| ![](https://i.imgur.com/vAGKjqx.png) |![](https://i.imgur.com/aWTMuUX.png) | ![](https://i.imgur.com/AL0vrco.png) |![](https://i.imgur.com/niEZNwX.png) | ![](https://i.imgur.com/HDDbiNz.png) |
| ![](https://i.imgur.com/tIQRU9A.png) |![](https://i.imgur.com/0FbXuB5.png) | ![](https://i.imgur.com/hwO6LVO.png) |![](https://i.imgur.com/gDinIVM.png) | ![](https://i.imgur.com/XBkavRr.png) |
| ![](https://i.imgur.com/zD0T7WV.png) |![](https://i.imgur.com/GHRpyiS.png) | ![](https://i.imgur.com/2Yl4VPg.png) |![](https://i.imgur.com/JvBPqVW.png) | ![](https://i.imgur.com/kA03Ukk.png) |
| ![](https://i.imgur.com/ACOkOTI.png) | ![](https://i.imgur.com/bbQNm4j.png) | ![](https://i.imgur.com/IxRVLG2.png) | ![](https://i.imgur.com/s804MLH.png) | ![](https://i.imgur.com/GrvdIEf.png) |



# [明龍的專題生們之專題紀錄](https://hackmd.io/drEMbOsRSwmAe4ulljQKFA?both)

## 側面 maximum intensity projection
![](https://i.imgur.com/ITvXzs9.png)
![](https://i.imgur.com/LKOKumq.png)
- 投影之後的結果如圖，在上面示意圖中的兩個病人的下排牙根以及下顎基本上已經很難分辨出來，經過幾次嘗試之後決定使用比較土法的方式(下限部分直接往下顎延伸50個 slices)
![](https://i.imgur.com/zj0mAmw.png =325x)![](https://i.imgur.com/JpFrmli.png =325x)

## Front MIP
1. 希望可以利用 Canny 找到牙齒的邊界以便切割有牙齒的部分。 :o: 
    - 如果希望牙根部分可以被識別的較為完整的話，圖片上會有很多的雜訊![](https://i.imgur.com/5oTThxc.png)
2. 參考[此論文](https://aip.scitation.org/doi/pdf/10.1063/1.4994425)，利用 Morphological Transformation 去除雜訊。
    - 去除雜訊後的 Canny 結果![](https://i.imgur.com/XkOeRPY.png)
    - 嘗試使用 GaussianBlur 發現結果更好一些![](https://i.imgur.com/KOGafWo.png)

3. 直接利用 threshold 二分法看能不能表示出牙齒 :x: 
    - 有些特別的 case 會難以處理，例如 case 8 整個 MIP 的兩度明顯比其他 case 高不少，因此沒有辦法順利找到牙齒的部分。![](https://i.imgur.com/WR5wtOQ.jpg)![](https://i.imgur.com/NWuhNbA.png)
    - 不過其實 Canny 在 case 8 的成效也不是很好![](https://i.imgur.com/HgPftjw.png)

- 目前使用方法：
    1. GaussianBlur
    2. Canny edge
    3. 計算 Canny output 的每個 row 的 255 數量，拿第一個以及最後一個超過 30 的中心點為基準，分別留下前後的 120 個 row。![](https://i.imgur.com/mPvMBEh.png)![](https://i.imgur.com/6dfCEt7.jpg)


## model
[以CNN切割單純分出牙齒](https://ieeexplore.ieee.org/document/9007457)
* metal artifact reduction (MAR)
    * sinograms that contain raw two-dimensional projections 
*  U-Net structure
* 由於資料量較少 -> multiphase learning method
* CBCT中牙齒體積很少 -> under-sample
* under-sample造成被裁切掉的部分容易出問題 -> multi-stage training
    1. 只含有牙齒的矩形(長寬高都做裁切)
    2. 有牙齒的slice(只切高)
    3. 完整CBCT影像
* 

[非機器學習分層](https://www.ncbi.nlm.nih.gov/pmc/articles/PMC6476381/)
* pre-processing stage做高斯轉換
* 以greyscale, neighbour average greyscale and gradient三種資訊來找出threshold 
    * using plane intercept histogram of reciprocal cross entropy algorithm

[以CNN切割單純分出牙齒](https://openaccess.thecvf.com/content_CVPR_2019/papers/Cui_ToothNet_Automatic_Tooth_Instance_Segmentation_and_Identification_From_Cone_Beam_CVPR_2019_paper.pdf)
[測試切割結果的準確度](http://oralpathol.dlearn.kmu.edu.tw/case/Journal%20reading-intern-17-5/CBCT-tooth%20structure-OOOR-2017.pdf)
[An Automatic Segmentation and Classification Framework Based on PCNN Model for Single Tooth in MicroCT Images](https://www.ncbi.nlm.nih.gov/pmc/articles/PMC4913904/)
[feature pyramid network](https://www.sciencedirect.com/science/article/pii/S0099239921000042?casa_token=T86X07zWWTkAAAAA:I4Ai40sqzSxgo7Sa_GJOdDUEYsQsEb2nZ7VldaQl2XE16ppnJds6Kfx5TXSWWssoKY5fXpLV6w)

---
GaussianMixtureModel Unsupervised![](https://i.imgur.com/rjskqQb.png)
1. 將整個影像的資料的 intensity 畫到用柱狀圖呈現如圖，經過測試之後確定 intensity 在 100 以下的皆為雜訊，所以我先使用 threshold 將雜訊去除。
2. 把非零的數據做標準化，使其平均分散在 0-255
3. 為了再加強資料的對比度，我還使用 gamma enhanced 增強資料的離散度
![](https://i.imgur.com/ky3KW2b.png)
4. 最後，使用 GMM 將資料分為三個類別

result:
![](https://i.imgur.com/WuugjgB.jpg)

![](https://i.imgur.com/aMJYabh.jpg)

---
## 使用UNet
* 程式碼來源説明：https://blog.csdn.net/u012931582/article/details/70215756
* 原理説明：https://blog.csdn.net/fjsd155/article/details/89107367
* github： https://github.com/zhixuhao/unet
* 使用的技術： 利用圖像扭曲技術，解決資料不足的問題
* input & output: png
* 目前的測試結果： 目前不清楚爲什麽出來的都是單色圖片（不一定是0，出現過極小數，刪去所有沒有label的圖片,可以發現到刪掉越多的圖片，最後的數字有提升，目前猜測有overfitting的問題）$\rightarrow$ 另有一種推測，是和loss function 有關，目前已經正在測試dice coeffient loss function. [各樣loss function 介紹](https://blog.csdn.net/m0_37477175/article/details/83004746)



## label 問題
牙齒中在看起來並非牙根的地方似乎有空洞，不知道醫生是只標記牙齒的外層，還是caseBai\DICOMExportTemp_20201202_172426\ScalarVolume_21中的檔案已經是模型預測出來的牙齒位置了
![](https://i.imgur.com/PN4FsxL.png)

# 2021/04/23
- 發現之前太早作標準化，導致 unsupervised 效果沒有很好，於是取消 normalization 並且換使用 hist equalization 使資料的分散度更明顯。
不過，目前的結果還是沒有辦法分辨下顎以及牙齒，目前的想法是不利用 processing 的方式是的數據更加分散，之後當作 unet feature 使用。![](https://i.imgur.com/9YNz2Ja.png)

- unsupervised 經過我的嘗試，讓 gmm 將資料分成 3 個 class 的結果比較好，不過也產生無法確定哪一個 class 是牙齒的 class。
（目前我的做法是：比較所有 class 的 max intensity，不過如果雜訊只要有一個點的亮度超過牙齒最高亮度得到的結果就會是錯誤的）

- unet 3D 結構完成，不過我家裡的電腦沒有足夠的記憶體，明天使用實驗室的電腦 run 看看目前已經有 lable 的 case

## U-Net目前結果
* 用訓練資料集的結果
    * * Original(IMG0184.png)
    * ![](https://i.imgur.com/OgTksye.png)
    * * Predict （Steps per Epoch = 2000. Epoch = 20）
    * ![](https://i.imgur.com/HuZpdHa.png)
    * Predict (Steps per Epoch = 2000, Epoch = 30)
    * ![](https://i.imgur.com/x5DmZ1y.png)

* 以case0035做測試（Steps per Epoch = 2000. Epoch = 20）
    * 原圖有牙齒的（效果比較好的）（170.png）
    * ![](https://i.imgur.com/5dOBgxa.png)
    * 預測的圖片
    * ![](https://i.imgur.com/31BN0JI.png)
    * 原圖有牙齒的（效果差的）(195.png)
    * ![](https://i.imgur.com/W9KosIv.png)
    * 預測的圖片
    * ![](https://i.imgur.com/paAF5RH.png)
    * 原圖無牙齒的(380.png)
    * ![](https://i.imgur.com/XIKQi34.png)
    * 預測的圖片
    * ![](https://i.imgur.com/Uob931c.png)

* 有發現一個狀況：如果是有植牙的情況，第一步需要先轉成位深8的png這一步，會出現不自然空洞的情況，目前還不清楚是什麽原因造成的
* 需要對閃光的部分做處理
* 無牙齒圖片的訓練/處理

# 2021/05/06
- model:
    - inputs: image(threshold) & canny edge
    - structure: unet with batch normalization
    - loss: binary_crossentropy
- data:
    - augmentation: flip & rotate
    - train data: 只是用完整的 case，且 input 為 dicom 檔

![](https://i.imgur.com/bUs9ouM.png)
![](https://i.imgur.com/cKTy5tW.png)


### post processing
* 離群的白點 -> 直接以座標來排除
* 磨平表面 -> 主要以morphology的opening來處理
![](https://i.imgur.com/MhQoLfl.png)
![](https://i.imgur.com/UDazB3m.png)
![](https://i.imgur.com/q4wbVhH.png)
![](https://i.imgur.com/kEjkYd0.png)
![](https://i.imgur.com/JVViajy.png)
![](https://i.imgur.com/cmmg8ah.png)

* 與牙齒接近的雜訊 -> 目前想以體積來處理
* 有看到一個以CRF來做後處理的，但由於需要再進行訓練，外加我看到別人再次實現的結果和論文差很多，所以先保留[Post Processing of Image Segmentation using Conditional Random Fields](https://ieeexplore.ieee.org/document/8991232)

# UNET
base_filter = 64

    inputs = Input(input_size)
    canny_input = Input(input_size)

    conv1 = Conv2D(base_filter, 3, padding='same', kernel_initializer='he_normal')(inputs)
    conv1 = Conv2D(base_filter, 3, padding='same', kernel_initializer='he_normal')(conv1)

    batch1 = BatchNormalization()(conv1)

    pool1 = MaxPooling2D(pool_size=(2, 2))(batch1)
    conv2 = Conv2D(base_filter*2, 3, padding='same', kernel_initializer='he_normal')(pool1)
    conv2 = Conv2D(base_filter*2, 3, padding='same', kernel_initializer='he_normal')(conv2)

    batch2 = BatchNormalization()(conv2)

    pool2 = MaxPooling2D(pool_size=(2, 2))(batch2)
    conv3 = Conv2D(base_filter*4, 3, padding='same', kernel_initializer='he_normal')(pool2)
    conv3 = Conv2D(base_filter*4, 3, padding='same', kernel_initializer='he_normal')(conv3)

    batch3 = BatchNormalization()(conv3)

    pool3 = MaxPooling2D(pool_size=(2, 2))(batch3)
    conv4 = Conv2D(base_filter*8, 3, padding='same', kernel_initializer='he_normal')(pool3)
    conv4 = Conv2D(base_filter*8, 3, padding='same', kernel_initializer='he_normal')(conv4)

    batch4 = BatchNormalization()(conv4)
    # drop4 = Dropout(0.5)(conv4)
    pool4 = MaxPooling2D(pool_size=(2, 2))(batch4)

    conv5 = Conv2D(base_filter*16, 3, padding='same', kernel_initializer='he_normal')(pool4)
    conv5 = Conv2D(base_filter*16, 3, padding='same', kernel_initializer='he_normal')(conv5)
    # drop5 = Dropout(0.5)(conv5)

    up6 = UpSampling2D(size=(2, 2))(conv5)
    merge6 = concatenate([conv4, up6], axis=3)
    conv6 = Conv2D(base_filter*8, 3, padding='same', kernel_initializer='he_normal')(merge6)
    conv6 = Conv2D(base_filter*8, 3, padding='same', kernel_initializer='he_normal')(conv6)

    up7 = UpSampling2D(size=(2, 2))(conv6)
    merge7 = concatenate([conv3, up7], axis=3)
    conv7 = Conv2D(base_filter*4, 3, padding='same', kernel_initializer='he_normal')(merge7)
    conv7 = Conv2D(base_filter*4, 3, padding='same', kernel_initializer='he_normal')(conv7)

    up8 = UpSampling2D(size=(2, 2))(conv7)
    merge8 = concatenate([conv2, up8], axis=3)
    conv8 = Conv2D(base_filter*2, 3, padding='same', kernel_initializer='he_normal')(merge8)
    conv8 = Conv2D(base_filter*2, 3, padding='same', kernel_initializer='he_normal')(conv8)

    up9 = UpSampling2D(size=(2, 2))(conv8)
    merge9 = concatenate([conv1, up9], axis=3)
    conv9 = Conv2D(base_filter, 3, padding='same', kernel_initializer='he_normal')(merge9)
    merge_last = concatenate([conv9, canny_input], axis=3)
    conv9 = Conv2D(base_filter, 3, padding='same', kernel_initializer='he_normal')(merge_last)
    conv10 = Conv2D(2, 3, padding='same', kernel_initializer='he_normal')(conv9)
    conv11 = Conv2D(1, 1, activation='sigmoid')(conv10)

    model = Model([inputs, canny_input], [conv11])

    # model.compile(optimizer='adam', loss='binary_crossentropy', metrics = [iou_acc])
    model.compile(optimizer=Adam(lr=1e-5), loss=dice_loss, metrics=['acc'])