PyTorch를 이용한 DICOM / EfficientNet B0 베이스라인 모델

ref_파이토치 + EfficientNetB0 Base + CV + 3 channel

Pytorch + EfficientNetB0 + CV 3 fold + 3 Channel

  • CV로 Submission 대신 정확도 측정
  • 3 channel로 정확도를 높임.
In [1]:
from apex import amp
import apex
import cv2
import glob
import pydicom
import numpy as np
import pandas as pd
from efficientnet_pytorch import EfficientNet
import torch
import torch.optim as optim
from albumentations import Compose, ShiftScaleRotate, Resize
from albumentations import CenterCrop, HorizontalFlip
from albumentations.pytorch import ToTensor
from torch.utils.data import Dataset
from tqdm import tqdm_notebook as tqdm
from matplotlib import pyplot as plt
import pydicom
import os
from sklearn.model_selection import StratifiedKFold
import datetime
from sklearn.metrics import log_loss
In [2]:
DESIRED_SIZE = (256,256)
N_CLASS = 6
N_EPOCH = 1
BATCH_SIZE = 32
N_SPLIT = 3

Dataset 구성

  • DICOM으로부터 바로 파일 구성

윈도우 스케일 교체

  • 스케일을 총 세개의 질병을 탐지할수 있게 바꿈
In [3]:
def win_scale(data, wl, ww, dtype, out_range):
    """
    Scale pixel intensity data using specified window level, width, and intensity range.
    """
    
    data_new = np.empty(data.shape, dtype=np.double)
    data_new.fill(out_range[1]-1)
    
    data_new[data <= (wl-ww/2.0)] = out_range[0]
    data_new[(data>(wl-ww/2.0))&(data<=(wl+ww/2.0))] = \
         ((data[(data>(wl-ww/2.0))&(data<=(wl+ww/2.0))]-(wl-0.5))/(ww-1.0)+0.5)*(out_range[1]-out_range[0])+out_range[0]
    data_new[data > (wl+ww/2.0)] = out_range[1]-1
    
    return data_new.astype(dtype)
        

def ct_win(im, wl, ww, dtype, out_range):    
    """
    Scale CT image represented as a `pydicom.dataset.FileDataset` instance.
    """

    # Convert pixel data from Houndsfield units to intensity:
    intercept = int(im[(0x0028, 0x1052)].value)
    slope = int(im[(0x0028, 0x1053)].value)
    data = (slope*im.pixel_array+intercept)

    # Scale intensity:
    return win_scale(data, wl, ww, dtype, out_range)
In [4]:
class HerDataset(Dataset):
    def __init__(self, target_data, path, labels, transform=None):
        self.path = path
        self.data = target_data
        self.transform = transform
        self.labels = labels
        
    def __len__(self):
        return len(self.data)
    
    def __getitem__(self, idx):
        
        img_name = os.path.join(self.path, self.data.loc[idx, 'Image'] + '.dcm')
        
#         img = self._read(img_name, DESIRED_SIZE)
        img = self._read_3c(img_name, DESIRED_SIZE)
        
        
        if self.transform:
            augment = self.transform(image = img)
            img = augmented['image']
            
        if self.labels:
            labels = torch.tensor(self.data.loc[idx, ['epidural', 'intraparenchymal', 'intraventricular', 'subarachnoid', 'subdural', 'any']])
            return {'image' : img, 'labels' : labels}
        else:
            return {'image' : img}
        
        
    def _read_3c(self, path, desired_size):
        
               
        ds = pydicom.dcmread(path)
        
        data1=ct_win(ds, 40,80, 'u1', (0,255))
        data2=ct_win(ds, 480,2500, 'u1', (0,255))
        data3=ct_win(ds, 90,350, 'u1', (0,255))
        
        data1 = cv2.resize(data1, desired_size, interpolation=cv2.INTER_LINEAR)
        data2 = cv2.resize(data2, desired_size, interpolation=cv2.INTER_LINEAR)
        data3 = cv2.resize(data3, desired_size, interpolation=cv2.INTER_LINEAR)
        
        return np.stack((data1, data2, data3))
        
        
    def _read(self, path, desired_size):
    
        dcm = pydicom.dcmread(path)

        window_params = self._get_windowing(dcm) # (center, width, slope, intercept)

        try:
            img = self._window_image(dcm.pixel_array, *window_params)
        except:
            img = np.zeros(desired_size)

        img = self._normalize(img)

        # 512, 512가 아니라면, 픽셀을 다시 맞춰줌 (적은 해상도, 큰 해상도 모두 Try할 필요 있음)
        
        img = cv2.resize(img, desired_size, interpolation=cv2.INTER_LINEAR)
    

        return np.squeeze(np.stack([img[np.newaxis,:,:]]*3, axis = 1))  #TODO, 3channel 때문에, 이걸 여러 디멘전으로 바꾼다   
    
    
    def _normalize(self, img):
        if img.max() == img.min():
            return np.zeros(img.shape)
        return 2 * (img - img.min())/(img.max() - img.min()) - 1
    
    def _window_image(self, img, window_center, window_width, slope, intercept):
        img = (img * slope + intercept)
        img_min = window_center - window_width//2
        img_max = window_center + window_width//2
        return np.clip(img, img_min, img_max) 
    
    def _get_windowing(self, data):
        dicom_fields = [data.WindowCenter, data.WindowWidth, data.RescaleSlope, data.RescaleIntercept]
        return [self._get_first_of_dicom_field_as_int(x) for x in dicom_fields]
    
    def _get_first_of_dicom_field_as_int(self, x):
        if type(x) == pydicom.multival.MultiValue:
            return int(x[0])
        else:
            return int(x)
    
    # TODO : save model epoch
        
    

로딩 및 세팅

  • EPOCH은 1~
  • N_SPLIT 은 CV나누는 기준
In [5]:
train = pd.read_csv('../input/stage_1_train.csv')
test = pd.read_csv('../input/stage_1_sample_submission.csv')

데이터 로딩

  • Augmentation 구성 / 여러가지를 돌려봐야함
In [6]:
transform_train = Compose([
#     ShiftScaleRotate(),
#     CenterCrop(200,200),
#     HorizontalFlip(),
    ToTensor()
])

transofrm_test = Compose([
#     CenterCrop(200,200),
#     HorizontalFlip(),
    ToTensor(),
])
  • 미리 구성된 pivot된 csv를 로딩한다
In [7]:
train_table = pd.read_csv('train_table.csv')
test_table = pd.read_csv('test.csv')

splits = list(StratifiedKFold(n_splits=N_SPLIT, shuffle=True, random_state=18).split(train_table['Image'], train_table['any']))


model_list = []
valid_preds = []
tr_losses = []

CV판단용 predict

  • CV에 대해서 log_loss를 계산하기 위한 predict 함수
In [8]:
def predict(model, dataset, dataloader):
    
    for param in model.parameters():
        param.requires_grad = False

    model.eval()

    test_pred = np.zeros((len(dataset) * N_CLASS, 1))

    for i, x_batch in enumerate(tqdm(dataloader)):

        x_batch = x_batch["image"]
        x_batch = x_batch.to(device, dtype=torch.float)

        with torch.no_grad():

            pred = model(x_batch)

            test_pred[(i * BATCH_SIZE * N_CLASS):((i + 1) * BATCH_SIZE * N_CLASS)] = torch.sigmoid(
                pred).detach().cpu().reshape((len(x_batch) * N_CLASS, 1))
    
    return test_pred

각각의 fold에 따른 CV구성

  • 3개로 잘라 작동.
In [ ]:
for sp in splits :
    train_target_table = train_table.iloc[sp[0]].copy().reset_index()
    valid_target_table = train_table.iloc[sp[1]].copy().reset_index()
    
    train_dataset = HerDataset(target_data = train_target_table, path = '../input/stage_1_train_images/', labels = True)
    valid_dataset = HerDataset(target_data = valid_target_table, path = '../input/stage_1_train_images/', labels = True)
    test_dataset = HerDataset(target_data = test_table, path = '../input/stage_1_test_images/', labels = False)

    # TODO : NUM_WORKER?
    data_loader_train = torch.utils.data.DataLoader(train_dataset, batch_size = BATCH_SIZE, shuffle = False, num_workers= 16)
    data_lader_valid = torch.utils.data.DataLoader(valid_dataset, batch_size = BATCH_SIZE, shuffle = False, num_workers= 16)
    data_loader_test = torch.utils.data.DataLoader(test_dataset, batch_size = BATCH_SIZE, shuffle = False, num_workers= 16)
    
    device = torch.device("cuda:0")
    model = EfficientNet.from_pretrained('efficientnet-b0')
    model._fc = torch.nn.Linear(1280, N_CLASS)
    model.to(device)

    criterion = torch.nn.BCEWithLogitsLoss()
    plist = [{'params' : model.parameters(), 'lr' : 2e-4}]
    optimizer = optim.Adam(plist, lr = 2e-4)

    model, optimizer = amp.initialize(model, optimizer, opt_level = "O1")
    
    for epoch in range(N_EPOCH):
        
        print('Epoch {} / {}'.format(epoch, N_EPOCH - 1))
        print('-' * 10)

        model.train()
        tr_loss = 0

        tk0 = tqdm(data_loader_train, desc="iteration")

        for step, batch in enumerate(tk0):
            inputs = batch['image']
            labels = batch['labels']

            inputs = inputs.to(device, dtype=torch.float)
            labels = labels.to(device, dtype=torch.float)

            outputs = model(inputs)
            loss = criterion(outputs, labels)

            with amp.scale_loss(loss, optimizer) as scaled_loss:
                scaled_loss.backward()

            tr_loss += loss.item()
            if (step % 100 == 0):
                print('traning loss : {:.4f}'.format(tr_loss / (step+1)))
                tr_losses.append(tr_loss / (step+1))
                
            optimizer.step()
            optimizer.zero_grad()

        valid_pred = predict(model, valid_dataset, data_lader_valid)    
        valid_preds.append(valid_pred)
            
        epoch_loss = tr_loss / len(data_loader_train)
        model_list.append(model)
        print('training loss : {:.4f}'.format(epoch_loss))
Loaded pretrained weights for efficientnet-b0
Selected optimization level O1:  Insert automatic casts around Pytorch functions and Tensor methods.

Defaults for this optimization level are:
enabled                : True
opt_level              : O1
cast_model_type        : None
patch_torch_functions  : True
keep_batchnorm_fp32    : None
master_weights         : None
loss_scale             : dynamic
Processing user overrides (additional kwargs that are not None)...
After processing overrides, optimization options are:
enabled                : True
opt_level              : O1
cast_model_type        : None
patch_torch_functions  : True
keep_batchnorm_fp32    : None
master_weights         : None
loss_scale             : dynamic
Warning:  multi_tensor_applier fused unscale kernel is unavailable, possibly because apex was installed without --cuda_ext --cpp_ext. Using Python fallback.  Original ImportError was: ModuleNotFoundError("No module named 'amp_C'")
Epoch 0 / 0
----------
traning loss : 0.7114
traning loss : 0.2491

각각의 CV 수렴 확인

  • 각각의 모델이 수렴하는지 확인
In [ ]:
# 모델을 하나만 돌리고 끊었을때
if len(model_list) == 1:

    valid_pred = predict(model, valid_dataset, data_lader_valid)    
    valid_preds.append(valid_pred)

    epoch_loss = tr_loss / len(data_loader_train)
    model_list.append(model)
    print('training loss : {:.4f}'.format(epoch_loss))
In [ ]:
import seaborn as sns; sns.set()
import matplotlib.pyplot as plt
ax = sns.lineplot(data=np.array(tr_losses))
In [ ]:
len(test_dataset)

Model 저장

  • model_list안에 총 세개의 모델이 있음
In [ ]:
import pickle

tm = str(datetime.datetime.now())[:16]

pickle_out = open("models_base" + tm + ".pickle","wb")
pickle.dump(model_list, pickle_out)
pickle_out.close()

Validation 정확도 측정

  • 3 cv로 구성. 시간이 굉장히 오래 걸려서 빠르게 측정하는 방법이 필요함.
In [ ]:
log_losses = []

for i in range(len(model_list)):
    
    val_log_losses = []
    
    for idx, loc in tqdm(enumerate(splits[i][1])):
        label = train_table.iloc[loc]
        current_idx = idx * 6        
        valid_perc = valid_preds[i][current_idx : current_idx + 6]
        
        
        ll = log_loss([label['epidural'], label['intraparenchymal'], label['intraventricular'], label['subarachnoid'], label['subdural'], label['any']]
            , valid_perc, sample_weight=([1,1,1,1,1,2]), labels=[0,1])
        
        val_log_losses.append(ll)
    
    print('valid_loss : ' + str(np.mean(val_log_losses)))
    
    log_losses = log_losses + val_log_losses
    
print('total_valid_loss : ' + str(np.mean(log_losses)))    

Test Prediction

  • 최종적인 결과는 총 3개의 결과를 각각 평균한다.
  • sigmoid 함수 안에서 평균치지 않도록 조심
In [ ]:
final_pred = np.zeros((len(test_dataset) * N_CLASS, 1))
sibal = []

for model in model_list:
    
    for param in model.parameters():
        param.requires_grad = False

    model.eval()

    test_pred = np.zeros((len(test_dataset) * N_CLASS, 1))

    for i, x_batch in enumerate(tqdm(data_loader_test)):

        x_batch = x_batch["image"]
        x_batch = x_batch.to(device, dtype=torch.float)

        with torch.no_grad():

            pred = model(x_batch)

            test_pred[(i * BATCH_SIZE * N_CLASS):((i + 1) * BATCH_SIZE * N_CLASS)] =  torch.sigmoid(
                pred).detach().cpu().reshape((len(x_batch) * N_CLASS, 1)) 
            
    sibal.append(test_pred)
            
    final_pred += ( test_pred / len(model_list) )
    
In [ ]:
len(final_pred)

Test 분포 확인

  • 일단 0.100의 제출용이라면, 평균은 0.05 근처, 표준편차는 0.15 근처로 보임.
In [ ]:
print(np.mean(final_pred))
print(np.std(final_pred))
print(len(final_pred))

제출용 파일 만들기

  • 시간에 따라 다르게 나뉨
  • any가 충분히 확률이 큰지 확인
In [ ]:
submission =  pd.read_csv('../input/stage_1_sample_submission.csv')
submission = pd.concat([submission.drop(columns=['Label']), pd.DataFrame(test_pred)], axis=1)
submission.columns = ['ID', 'Label']

print(str(tm))

submission.to_csv('submission_b0' + tm + '.csv', index=False)
submission.head(10)

답글 남기기