{"cells":[{"metadata":{},"cell_type":"markdown","source":"# OSIC: keras model inference on images and tabular data","execution_count":null},{"metadata":{},"cell_type":"markdown","source":"This is an inference part for a notebook with model training [OSIC keras images and tabular data model](https://www.kaggle.com/vgarshin/osic-keras-images-and-tabular-data-model).","execution_count":null},{"metadata":{"trusted":true},"cell_type":"code","source":"%%time\n!pip install ../input/kerasapplications/keras-team-keras-applications-3b180cb -f ./ --no-index\n!pip install ../input/efficientnet/efficientnet-1.1.0/ -f ./ --no-index","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","trusted":true},"cell_type":"code","source":"import warnings\nwarnings.simplefilter(action='ignore', category=FutureWarning)\nimport numpy as np\nimport pandas as pd\nimport os\nimport cv2\nimport random\nimport matplotlib.pyplot as plt\nfrom tqdm.notebook import tqdm\nfrom sklearn.preprocessing import MinMaxScaler\nfrom sklearn.model_selection import KFold, GroupKFold\nfrom skimage.segmentation import clear_border\nfrom skimage.measure import label, regionprops\nimport tensorflow as tf\nimport tensorflow.keras.backend as K\nimport tensorflow.keras.layers as L\nimport tensorflow.keras.models as M\nimport efficientnet.tfkeras as efn\nfrom tensorflow.keras.utils import Sequence\nfrom tensorflow.keras.optimizers import Adam\nfrom tensorflow.keras.callbacks import EarlyStopping, LearningRateScheduler, ModelCheckpoint, ReduceLROnPlateau\nfrom tensorflow.keras.models import load_model\nimport pydicom\nprint('tensorflow version:', tf.__version__)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def seed_all(seed=2020):\n    random.seed(seed)\n    os.environ['PYTHONHASHSEED'] = str(seed)\n    np.random.seed(seed)\n    tf.random.set_seed(seed)\nseed_all(2020)\nDATA_PATH = '../input/osic-pulmonary-fibrosis-progression'\nBATCH_SIZE_PRED = 1\nFEATURES = True\nADV_FEATURES = False\nC_SIGMA, C_DELTA = tf.constant(70, dtype='float32'), tf.constant(1000, dtype='float32')\nQS = [.05, .50, .95]\nIMG_SIZE = 224\nRESIZE = 224\nSEQ_LEN = 12\nCUTOFF = 2\nLAMBDA = .8\nMDL_VERSION = 'v5'\nMODELS_PATH = '../input/osic-keras-trained-model-v0'\nMASKED = False","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train = pd.read_csv(f'{DATA_PATH}/train.csv')\ntrain.drop_duplicates(keep=False, inplace=True, subset=['Patient','Weeks'])\ntest = pd.read_csv(f'{DATA_PATH}/test.csv')\nsubm = pd.read_csv(f'{DATA_PATH}/sample_submission.csv')\nsubm['Patient'] = subm['Patient_Week'].apply(lambda x: x.split('_')[0])\nsubm['Weeks'] = subm['Patient_Week'].apply(lambda x: int(x.split('_')[-1]))\nsubm =  subm[['Patient','Weeks','Confidence','Patient_Week']]\nsubm = subm.merge(test.drop('Weeks', axis=1), on='Patient')\ntrain['SPLIT'] = 'train'\ntest['SPLIT'] = 'val'\nsubm['SPLIT'] = 'test'\ndata = train.append([test, subm])\nprint('train:',  train.shape, 'unique Pats:', train.Patient.nunique(),\n      '\\ntest:', test.shape,  'unique Pats:', test.Patient.nunique(),\n      '\\nsubm:', subm.shape,  'unique Pats:', subm.Patient.nunique(),\n      '\\ndata',  data.shape,  'unique Pats:', data.Patient.nunique())\ndata['min_week'] = data['Weeks']\ndata.loc[data.SPLIT == 'test', 'min_week'] = np.nan\ndata['min_week'] = data.groupby('Patient')['min_week'].transform('min')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"data = pd.concat([data, pd.get_dummies(data.Sex), pd.get_dummies(data.SmokingStatus)], axis=1)\nif FEATURES:\n    base = data.loc[data.Weeks == data.min_week]\n    base = base[['Patient', 'FVC']].copy()\n    base.columns = ['Patient', 'min_week_FVC']\n    base['nb'] = 1\n    base['nb'] = base.groupby('Patient')['nb'].transform('cumsum')\n    base = base[base.nb == 1]\n    base.drop('nb', axis=1, inplace=True)\n    data = data.merge(base, on='Patient', how='left')\n    data['relative_week'] = data['Weeks'] - data['min_week']\n    del base\nif ADV_FEATURES:\n    target_cols = ['FVC']\n    enc_cols =  [ \n        'Female',\n        'Male',\n        'Currently smokes',\n        'Ex-smoker',\n        'Never smoked'\n    ]\n    for t_col in target_cols:\n        for col in enc_cols:\n            col_name = f'_{col}_{t_col}_'\n            data[f'enc{col_name}mean'] = data.groupby(col)[t_col].transform('mean')\n            data[f'enc{col_name}std'] = data.groupby(col)[t_col].transform('std')\n    data['TC'] = 0\n    data.loc[data['Weeks'] == 0, 'TC'] = 1\nprint(data.shape)\nprint(data.columns)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"feat_cols = [\n    'Female', 'Male',\n    'Currently smokes', \n    'Ex-smoker', 'Never smoked'\n]\nscale_cols = [\n    'Percent', \n    'Age', \n    'relative_week', \n    'min_week_FVC'\n]\nscale_cols.extend([x for x in data.columns if 'FVC_mean' in x])\nscale_cols.extend([x for x in data.columns if 'FVC_std' in x])\nscaler = MinMaxScaler()\ndata[scale_cols] = scaler.fit_transform(data[scale_cols])\nfeat_cols.extend(scale_cols)\nprint('all data columns:', data.columns)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train = data.loc[data.SPLIT == 'train']\ntest = data.loc[data.SPLIT == 'val']\nsubm = data.loc[data.SPLIT == 'test']\ndel data\nsubm.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"class DataGenOsic(Sequence):\n    def __init__(self, df, tab_cols,\n                 batch_size=8, mode='fit', shuffle=False, \n                 aug=None, resize=None, masked=True, cutoff=2,\n                 seq_len=12, img_size=224):\n        self.df = df\n        self.shuffle = shuffle\n        self.mode = mode\n        self.aug = aug\n        self.resize = resize\n        self.masked = masked\n        self.cutoff = cutoff\n        self.batch_size = batch_size\n        self.img_size = img_size\n        self.seq_len = seq_len\n        self.tab_cols = tab_cols\n        self.on_epoch_end()\n    def __len__(self):\n        return int(np.floor(len(self.df) / self.batch_size))\n    def on_epoch_end(self):\n        self.indexes = np.arange(len(self.df))\n        if self.shuffle:\n            np.random.shuffle(self.indexes)\n    def __getitem__(self, index):\n        batch_size = min(self.batch_size, len(self.df) - index * self.batch_size)\n        X_img = np.zeros((batch_size, self.seq_len, self.img_size, self.img_size, 3), dtype=np.float32)\n        X_tab = self.df[index * self.batch_size : (index + 1) * self.batch_size][self.tab_cols].values\n        pats_batch = self.df[index * self.batch_size : (index + 1) * self.batch_size]['Patient'].values\n        for i, pat_id in enumerate(pats_batch):\n            imgs_seq = self.get_imgs_seq(pat_id)\n            X_img[i, ] = imgs_seq\n        if self.mode == 'fit':\n            y = np.array(\n                self.df[index * self.batch_size : (index + 1) * self.batch_size]['FVC'].values, \n                dtype=np.float32\n            )\n            return (X_img, X_tab), y\n        elif self.mode == 'predict':\n            y = np.zeros(batch_size, dtype=np.float32)\n            return (X_img, X_tab), y\n        else:\n            raise AttributeError('mode parameter error')\n    def load_scan(self, pat_id):\n        if self.mode == 'fit':\n            path = f'{DATA_PATH}/train/{pat_id}'\n        elif self.mode == 'predict':\n            path = f'{DATA_PATH}/test/{pat_id}'\n        else:\n            raise AttributeError('mode parameter error')\n        file_names = sorted(os.listdir(path), key=lambda x: int(os.path.splitext(x)[0]))\n        idxs = [\n            int(i * len(file_names) / (self.seq_len + 2 * self.cutoff)) \n            for i in range(self.seq_len + 2 * self.cutoff)\n        ]\n        slices = [\n            pydicom.read_file(path + '/' + file_names[idx])\n            for idx in idxs[self.cutoff:-self.cutoff]\n        ]\n        if len(slices) < self.seq_len:\n            for i in range(self.seq_len - len(slices)):\n                slices.append(\n                    pydicom.read_file(path + '/' + os.listdir(path)[-1])\n                )\n        return slices\n    def get_pixels_hu(self, scans):\n        image = np.stack([s.pixel_array.astype(float) for s in scans])\n        image = image.astype(np.int16)\n        image[image == -2000] = 0\n        intercept = scans[0].RescaleIntercept\n        slope = scans[0].RescaleSlope\n        if slope != 1:\n            image = slope * image.astype(np.float64)\n            image = image.astype(np.int16)\n        image += np.int16(intercept)\n        return np.array(image, dtype=np.int16)\n    def get_imgs_seq(self, pat_id):\n        seq_imgs = []\n        slices = self.load_scan(pat_id)\n        scans = self.get_pixels_hu(slices)\n        for img_idx in range(self.seq_len):\n            img = scans[img_idx]\n            if self.masked:\n                mask = self.get_lungs_mask(img)\n                img[mask == False] = img.min()\n            if self.resize:\n                img = cv2.resize(img, (self.resize, self.resize))\n            img = (img - np.min(img)) / (np.max(img) - np.min(img))\n            img = np.repeat(img[..., np.newaxis], 3, -1)\n            seq_imgs.append(img)                 \n        return np.array(seq_imgs).astype(np.float32)\n    def get_lungs_mask(self, img):\n        mask = img < -200\n        mask = clear_border(mask)\n        mask = label(mask)\n        areas = [r.area for r in regionprops(mask)]\n        areas.sort()\n        if len(areas) > 2:\n            for region in regionprops(mask):\n                if region.area < areas[-2]:\n                    for coordinates in region.coords:\n                        mask[coordinates[0], coordinates[1]] = 0\n        return mask > 0","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def metric(y_true, y_pred, pred_std):\n    clip_std = np.clip(pred_std, 70, 9e9)  \n    delta = np.clip(np.abs(y_true - y_pred), 0 , 1000)  \n    return np.mean(-1 * (np.sqrt(2) * delta / clip_std) - np.log(np.sqrt(2) * clip_std))\ndef score(y_true, y_pred):\n    tf.dtypes.cast(y_true, tf.float32)\n    tf.dtypes.cast(y_pred, tf.float32)\n    sigma = (y_pred[:, 2] - y_pred[:, 0]) / 2\n    fvc_pred = y_pred[:, 1]\n    sigma_clip = tf.maximum(sigma, C_SIGMA)\n    delta = tf.abs(y_true[:, 0] - fvc_pred)\n    delta = tf.minimum(delta, C_DELTA)\n    sq2 = tf.sqrt(tf.dtypes.cast(2, dtype=tf.float32))\n    metric = sq2 * (delta / sigma_clip) + tf.math.log(sigma_clip * sq2)\n    return K.mean(metric)\ndef qloss(y_true, y_pred):\n    q = tf.constant(np.array([QS]), dtype=tf.float32)\n    e = y_true - y_pred\n    v = tf.maximum(q * e, (q - 1) * e)\n    return K.mean(v)\ndef mloss(lmbd):\n    def loss(y_true, y_pred):\n        return lmbd * qloss(y_true, y_pred) + (1 - lmbd) * score(y_true, y_pred)\n    return loss","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"%%time\nmodel_file = f'{MODELS_PATH}/model_{MDL_VERSION}.h5'\nmodel = load_model(model_file, custom_objects={'qloss': qloss, 'loss': mloss(LAMBDA), 'score': score})\nprint('model loaded:', model_file)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"**NOTE:** Parameter `masked=True` slows down prediction generator and causes `Notebook timeout` error. ","execution_count":null},{"metadata":{"trusted":true},"cell_type":"code","source":"subm_datagen = DataGenOsic(\n    df=subm,\n    tab_cols=feat_cols,\n    batch_size=BATCH_SIZE_PRED,\n    mode='predict', \n    shuffle=False, \n    aug=None,\n    resize=RESIZE,\n    masked=MASKED,\n    cutoff=CUTOFF,\n    seq_len=SEQ_LEN, \n    img_size=IMG_SIZE\n)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"(Xt_img, Xt_tab), _ = subm_datagen.__getitem__(0)\nprint('val X img: ', Xt_img.shape)\nprint('val X tab: ', Xt_tab.shape)\nfig, axes = plt.subplots(figsize=(10, 4), nrows=BATCH_SIZE_PRED, ncols=SEQ_LEN)\nfor j in range(BATCH_SIZE_PRED):\n    for i in range(SEQ_LEN):\n        if BATCH_SIZE_PRED > 1:\n            axes[j, i].imshow(Xt_img[j][i])\n            axes[j, i].axis('off')\n            axes[j, i].set_title(f'{j + 1} of {BATCH_SIZE_PRED}')\n        else:\n            axes[i].imshow(Xt_img[j][i])\n            axes[i].axis('off')\n            axes[i].set_title(f'{j + 1} of {BATCH_SIZE_PRED}')\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"%%time\npreds_subm = model.predict_generator(subm_datagen, verbose=1)\nprint('predictions shape:', preds_subm.shape)\nprint('predictions sample:', preds_subm[0])","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"subm['FVC'] = preds_subm[:, 1]\nsubm['Confidence'] = (preds_subm[:, 2] - preds_subm[:, 0]) / 2\nsubm[['Patient_Week','FVC','Confidence']].to_csv('submission.csv', index=False)\nsubm[['Patient_Week','FVC','Confidence']].describe().T","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"","execution_count":null,"outputs":[]}],"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":4,"nbformat_minor":4}