{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"import os\nimport cv2\nimport glob\nimport PIL\nimport shutil\nimport numpy as np\nimport pandas as pd\nimport seaborn as sns\nimport matplotlib.pyplot as plt\nfrom skimage import data\nfrom skimage.util import montage \nimport skimage.transform as skTrans\nfrom skimage.transform import rotate\nfrom skimage.transform import resize\nfrom PIL import Image, ImageOps  \n\n\n# neural imaging\nimport nilearn as nl\nimport nibabel as nib\nimport nilearn.plotting as nlplt\n!pip install git+https://github.com/miykael/gif_your_nifti # nifti to gif \nimport gif_your_nifti.core as gif2nif\n\n\n# ml libs\nimport keras\nimport keras.backend as K\nfrom keras.callbacks import CSVLogger\nimport tensorflow as tf\nfrom tensorflow.keras.utils import plot_model\nfrom sklearn.preprocessing import MinMaxScaler\nfrom sklearn.model_selection import train_test_split\nfrom sklearn.metrics import classification_report\nfrom tensorflow.keras.models import *\nfrom tensorflow.keras.layers import *\nfrom tensorflow.keras.optimizers import *\nfrom tensorflow.keras.callbacks import ModelCheckpoint, ReduceLROnPlateau, EarlyStopping, TensorBoard\nfrom tensorflow.keras.layers.experimental import preprocessing\n\n\n# Make numpy printouts easier to read.\nnp.set_printoptions(precision=3, suppress=True)","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","_kg_hide-output":true,"execution":{"iopub.status.busy":"2021-12-12T09:08:09.754290Z","iopub.execute_input":"2021-12-12T09:08:09.754742Z","iopub.status.idle":"2021-12-12T09:08:34.725086Z","shell.execute_reply.started":"2021-12-12T09:08:09.754653Z","shell.execute_reply":"2021-12-12T09:08:34.724013Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# DEFINE seg-areas  \nSEGMENT_CLASSES = {\n    0 : 'NOT tumor',\n    1 : 'NECROTIC/CORE', # or NON-ENHANCING tumor CORE\n    2 : 'EDEMA',\n    3 : 'ENHANCING' # original 4 -> converted into 3 later\n}\n\n# there are 155 slices per volume\n# to start at 5 and use 145 slices means we will skip the first 5 and last 5 \nVOLUME_SLICES = 100 \nVOLUME_START_AT = 22 # first slice of volume that we will include","metadata":{"execution":{"iopub.status.busy":"2021-12-12T09:08:34.727669Z","iopub.execute_input":"2021-12-12T09:08:34.728107Z","iopub.status.idle":"2021-12-12T09:08:34.733451Z","shell.execute_reply.started":"2021-12-12T09:08:34.728060Z","shell.execute_reply":"2021-12-12T09:08:34.731988Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"TRAIN_DATASET_PATH = '../input/brats20-dataset-training-validation/BraTS2020_TrainingData/MICCAI_BraTS2020_TrainingData/'\nVALIDATION_DATASET_PATH = '../input/brats20-dataset-training-validation/BraTS2020_ValidationData/MICCAI_BraTS2020_ValidationData'\n\ntest_image_flair=nib.load(TRAIN_DATASET_PATH + 'BraTS20_Training_001/BraTS20_Training_001_flair.nii').get_fdata()\ntest_image_t1=nib.load(TRAIN_DATASET_PATH + 'BraTS20_Training_001/BraTS20_Training_001_t1.nii').get_fdata()\ntest_image_t1ce=nib.load(TRAIN_DATASET_PATH + 'BraTS20_Training_001/BraTS20_Training_001_t1ce.nii').get_fdata()\ntest_image_t2=nib.load(TRAIN_DATASET_PATH + 'BraTS20_Training_001/BraTS20_Training_001_t2.nii').get_fdata()\ntest_mask=nib.load(TRAIN_DATASET_PATH + 'BraTS20_Training_001/BraTS20_Training_001_seg.nii').get_fdata()\n\n\nfig, (ax1, ax2, ax3, ax4, ax5) = plt.subplots(1,5, figsize = (20, 10))\nslice_w = 25\nax1.imshow(test_image_flair[:,:,test_image_flair.shape[0]//2-slice_w], cmap = 'gray')\nax1.set_title('Image flair')\nax2.imshow(test_image_t1[:,:,test_image_t1.shape[0]//2-slice_w], cmap = 'gray')\nax2.set_title('Image t1')\nax3.imshow(test_image_t1ce[:,:,test_image_t1ce.shape[0]//2-slice_w], cmap = 'gray')\nax3.set_title('Image t1ce')\nax4.imshow(test_image_t2[:,:,test_image_t2.shape[0]//2-slice_w], cmap = 'gray')\nax4.set_title('Image t2')\nax5.imshow(test_mask[:,:,test_mask.shape[0]//2-slice_w])\nax5.set_title('Mask')\n","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2021-12-12T09:08:34.734712Z","iopub.execute_input":"2021-12-12T09:08:34.735273Z","iopub.status.idle":"2021-12-12T09:08:36.443647Z","shell.execute_reply.started":"2021-12-12T09:08:34.735237Z","shell.execute_reply":"2021-12-12T09:08:36.442380Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# dice loss as defined above for 4 classes\ndef dice_coef(y_true, y_pred, smooth=1.0):\n    class_num = 4\n    for i in range(class_num):\n        y_true_f = K.flatten(y_true[:,:,:,i])\n        y_pred_f = K.flatten(y_pred[:,:,:,i])\n        intersection = K.sum(y_true_f * y_pred_f)\n        loss = ((2. * intersection + smooth) / (K.sum(y_true_f) + K.sum(y_pred_f) + smooth))\n   #     K.print_tensor(loss, message='loss value for class {} : '.format(SEGMENT_CLASSES[i]))\n        if i == 0:\n            total_loss = loss\n        else:\n            total_loss = total_loss + loss\n    total_loss = total_loss / class_num\n#    K.print_tensor(total_loss, message=' total dice coef: ')\n    return total_loss\n\n\n \n# define per class evaluation of dice coef\n# inspired by https://github.com/keras-team/keras/issues/9395\ndef dice_coef_necrotic(y_true, y_pred, epsilon=1e-6):\n    intersection = K.sum(K.abs(y_true[:,:,:,1] * y_pred[:,:,:,1]))\n    return (2. * intersection) / (K.sum(K.square(y_true[:,:,:,1])) + K.sum(K.square(y_pred[:,:,:,1])) + epsilon)\n\ndef dice_coef_edema(y_true, y_pred, epsilon=1e-6):\n    intersection = K.sum(K.abs(y_true[:,:,:,2] * y_pred[:,:,:,2]))\n    return (2. * intersection) / (K.sum(K.square(y_true[:,:,:,2])) + K.sum(K.square(y_pred[:,:,:,2])) + epsilon)\n\ndef dice_coef_enhancing(y_true, y_pred, epsilon=1e-6):\n    intersection = K.sum(K.abs(y_true[:,:,:,3] * y_pred[:,:,:,3]))\n    return (2. * intersection) / (K.sum(K.square(y_true[:,:,:,3])) + K.sum(K.square(y_pred[:,:,:,3])) + epsilon)\n\n\n\n# Computing Precision \ndef precision(y_true, y_pred):\n        true_positives = K.sum(K.round(K.clip(y_true * y_pred, 0, 1)))\n        predicted_positives = K.sum(K.round(K.clip(y_pred, 0, 1)))\n        precision = true_positives / (predicted_positives + K.epsilon())\n        return precision\n\n    \n# Computing Sensitivity      \ndef sensitivity(y_true, y_pred):\n    true_positives = K.sum(K.round(K.clip(y_true * y_pred, 0, 1)))\n    possible_positives = K.sum(K.round(K.clip(y_true, 0, 1)))\n    return true_positives / (possible_positives + K.epsilon())\n\n\n# Computing Specificity\ndef specificity(y_true, y_pred):\n    true_negatives = K.sum(K.round(K.clip((1-y_true) * (1-y_pred), 0, 1)))\n    possible_negatives = K.sum(K.round(K.clip(1-y_true, 0, 1)))\n    return true_negatives / (possible_negatives + K.epsilon())","metadata":{"execution":{"iopub.status.busy":"2021-12-12T09:08:36.445326Z","iopub.execute_input":"2021-12-12T09:08:36.445664Z","iopub.status.idle":"2021-12-12T09:08:36.465397Z","shell.execute_reply.started":"2021-12-12T09:08:36.445625Z","shell.execute_reply":"2021-12-12T09:08:36.464481Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"IMG_SIZE=128","metadata":{"execution":{"iopub.status.busy":"2021-12-12T09:08:36.467035Z","iopub.execute_input":"2021-12-12T09:08:36.467366Z","iopub.status.idle":"2021-12-12T09:08:36.482059Z","shell.execute_reply.started":"2021-12-12T09:08:36.467338Z","shell.execute_reply":"2021-12-12T09:08:36.480654Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# source https://naomi-fridman.medium.com/multi-class-image-segmentation-a5cc671e647a\n\ndef build_unet(inputs, ker_init, dropout):\n    conv1 = Conv2D(32, 3, activation = 'relu', padding = 'same', kernel_initializer = ker_init)(inputs)\n    conv1 = Conv2D(32, 3, activation = 'relu', padding = 'same', kernel_initializer = ker_init)(conv1)\n    \n    pool = MaxPooling2D(pool_size=(2, 2))(conv1)\n    conv = Conv2D(64, 3, activation = 'relu', padding = 'same', kernel_initializer = ker_init)(pool)\n    conv = Conv2D(64, 3, activation = 'relu', padding = 'same', kernel_initializer = ker_init)(conv)\n    \n    pool1 = MaxPooling2D(pool_size=(2, 2))(conv)\n    conv2 = Conv2D(128, 3, activation = 'relu', padding = 'same', kernel_initializer = ker_init)(pool1)\n    conv2 = Conv2D(128, 3, activation = 'relu', padding = 'same', kernel_initializer = ker_init)(conv2)\n    \n    pool2 = MaxPooling2D(pool_size=(2, 2))(conv2)\n    conv3 = Conv2D(256, 3, activation = 'relu', padding = 'same', kernel_initializer = ker_init)(pool2)\n    conv3 = Conv2D(256, 3, activation = 'relu', padding = 'same', kernel_initializer = ker_init)(conv3)\n    \n    \n    pool4 = MaxPooling2D(pool_size=(2, 2))(conv3)\n    conv5 = Conv2D(512, 3, activation = 'relu', padding = 'same', kernel_initializer = ker_init)(pool4)\n    conv5 = Conv2D(512, 3, activation = 'relu', padding = 'same', kernel_initializer = ker_init)(conv5)\n    drop5 = Dropout(dropout)(conv5)\n\n    up7 = Conv2D(256, 2, activation = 'relu', padding = 'same', kernel_initializer = ker_init)(UpSampling2D(size = (2,2))(drop5))\n    merge7 = concatenate([conv3,up7], axis = 3)\n    conv7 = Conv2D(256, 3, activation = 'relu', padding = 'same', kernel_initializer = ker_init)(merge7)\n    conv7 = Conv2D(256, 3, activation = 'relu', padding = 'same', kernel_initializer = ker_init)(conv7)\n\n    up8 = Conv2D(128, 2, activation = 'relu', padding = 'same', kernel_initializer = ker_init)(UpSampling2D(size = (2,2))(conv7))\n    merge8 = concatenate([conv2,up8], axis = 3)\n    conv8 = Conv2D(128, 3, activation = 'relu', padding = 'same', kernel_initializer = ker_init)(merge8)\n    conv8 = Conv2D(128, 3, activation = 'relu', padding = 'same', kernel_initializer = ker_init)(conv8)\n\n    up9 = Conv2D(64, 2, activation = 'relu', padding = 'same', kernel_initializer = ker_init)(UpSampling2D(size = (2,2))(conv8))\n    merge9 = concatenate([conv,up9], axis = 3)\n    conv9 = Conv2D(64, 3, activation = 'relu', padding = 'same', kernel_initializer = ker_init)(merge9)\n    conv9 = Conv2D(64, 3, activation = 'relu', padding = 'same', kernel_initializer = ker_init)(conv9)\n    \n    up = Conv2D(32, 2, activation = 'relu', padding = 'same', kernel_initializer = ker_init)(UpSampling2D(size = (2,2))(conv9))\n    merge = concatenate([conv1,up], axis = 3)\n    conv = Conv2D(32, 3, activation = 'relu', padding = 'same', kernel_initializer = ker_init)(merge)\n    conv = Conv2D(32, 3, activation = 'relu', padding = 'same', kernel_initializer = ker_init)(conv)\n    \n    conv10 = Conv2D(4, (1,1), activation = 'softmax')(conv)\n    \n    return Model(inputs = inputs, outputs = conv10)\n\ninput_layer = Input((IMG_SIZE, IMG_SIZE, 2))\n\nmodel = build_unet(input_layer, 'he_normal', 0.2)\nmodel.compile(loss=\"categorical_crossentropy\", optimizer=keras.optimizers.Adam(learning_rate=0.001))\nmodel.load_weights('../input/jkkghgyu/model_x1_1.h5')","metadata":{"_kg_hide-input":true,"_kg_hide-output":true,"execution":{"iopub.status.busy":"2021-12-12T09:08:43.370749Z","iopub.execute_input":"2021-12-12T09:08:43.371148Z","iopub.status.idle":"2021-12-12T09:08:44.896170Z","shell.execute_reply.started":"2021-12-12T09:08:43.371110Z","shell.execute_reply":"2021-12-12T09:08:44.895028Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"class DataGenerator(keras.utils.Sequence):\n    'Generates data for Keras'\n    def __init__(self, list_IDs, dim=(IMG_SIZE,IMG_SIZE), batch_size = 1, n_channels = 2, shuffle=True):\n        'Initialization'\n        self.dim = dim\n        self.batch_size = batch_size\n        self.list_IDs = list_IDs\n        self.n_channels = n_channels\n        self.shuffle = shuffle\n        self.on_epoch_end()\n\n    def __len__(self):\n        'Denotes the number of batches per epoch'\n        return int(np.floor(len(self.list_IDs) / self.batch_size))\n\n    def __getitem__(self, index):\n        'Generate one batch of data'\n        # Generate indexes of the batch\n        indexes = self.indexes[index*self.batch_size:(index+1)*self.batch_size]\n\n        # Find list of IDs\n        Batch_ids = [self.list_IDs[k] for k in indexes]\n\n        # Generate data\n        X, y = self.__data_generation(Batch_ids)\n\n        return X, y\n\n    def on_epoch_end(self):\n        'Updates indexes after each epoch'\n        self.indexes = np.arange(len(self.list_IDs))\n        if self.shuffle == True:\n            np.random.shuffle(self.indexes)\n\n    def __data_generation(self, Batch_ids):\n        'Generates data containing batch_size samples' # X : (n_samples, *dim, n_channels)\n        # Initialization\n        X = np.zeros((self.batch_size*VOLUME_SLICES, *self.dim, self.n_channels))\n        y = np.zeros((self.batch_size*VOLUME_SLICES, 240, 240))\n        Y = np.zeros((self.batch_size*VOLUME_SLICES, *self.dim, 4))\n\n        \n        # Generate data\n        for c, i in enumerate(Batch_ids):\n            case_path = os.path.join(TRAIN_DATASET_PATH, i)\n\n            data_path = os.path.join(case_path, f'{i}_flair.nii');\n            flair = nib.load(data_path).get_fdata()    \n\n            data_path = os.path.join(case_path, f'{i}_t1ce.nii');\n            ce = nib.load(data_path).get_fdata()\n            \n            data_path = os.path.join(case_path, f'{i}_seg.nii');\n            seg = nib.load(data_path).get_fdata()\n        \n            for j in range(VOLUME_SLICES):\n                 X[j +VOLUME_SLICES*c,:,:,0] = cv2.resize(flair[:,:,j+VOLUME_START_AT], (IMG_SIZE, IMG_SIZE));\n                 X[j +VOLUME_SLICES*c,:,:,1] = cv2.resize(ce[:,:,j+VOLUME_START_AT], (IMG_SIZE, IMG_SIZE));\n\n                 y[j +VOLUME_SLICES*c] = seg[:,:,j+VOLUME_START_AT];\n                    \n        # Generate masks\n        y[y==4] = 3;\n        mask = tf.one_hot(y, 4);\n        Y = tf.image.resize(mask, (IMG_SIZE, IMG_SIZE));\n        return X/np.max(X), Y\n        ","metadata":{"_kg_hide-input":false,"_kg_hide-output":true,"execution":{"iopub.status.busy":"2021-12-12T09:08:44.898670Z","iopub.execute_input":"2021-12-12T09:08:44.899129Z","iopub.status.idle":"2021-12-12T09:08:44.918911Z","shell.execute_reply.started":"2021-12-12T09:08:44.899082Z","shell.execute_reply":"2021-12-12T09:08:44.917810Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from tqdm import tqdm\ntest=DataGenerator(np.load('../input/hyfyuf/test (2).npy'))\ndef dice(t,p):\n    p=p/np.max(p)\n    p[p>0.8]=1\n    p[p<0.8]=0\n    return 2*np.sum(p*t)/(np.sum(t)+np.sum(p))\ndk={1:[],2:[],3:[]}\nfor r in tqdm(range(len(test))):\n    a=test[r]\n    p=model(a[0]).numpy()\n    for i in range(100):\n        for j in range(1,4):\n            true=a[1][i,:,:,j].numpy()\n            predicted=p[i,:,:,j]\n            if np.sum(true)>0:\n                dk[j].append(dice(true.copy(),predicted.copy()))","metadata":{"execution":{"iopub.status.busy":"2021-12-12T09:29:40.358380Z","iopub.execute_input":"2021-12-12T09:29:40.358931Z","iopub.status.idle":"2021-12-12T09:36:54.132176Z","shell.execute_reply.started":"2021-12-12T09:29:40.358898Z","shell.execute_reply":"2021-12-12T09:36:54.131242Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"np.mean(dk[1]),np.mean(dk[2]),np.mean(dk[3])","metadata":{"execution":{"iopub.status.busy":"2021-12-12T09:37:28.851158Z","iopub.execute_input":"2021-12-12T09:37:28.851866Z","iopub.status.idle":"2021-12-12T09:37:28.861474Z","shell.execute_reply.started":"2021-12-12T09:37:28.851813Z","shell.execute_reply":"2021-12-12T09:37:28.860429Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}