{"cells":[{"metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true},"cell_type":"code","source":"# %% [markdown]\n# # OSIC: keras model inference on images and tabular data\n\n# %% [markdown]\n# 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).\n\n# %% [code]\n# %%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\n\n# %% [code]\nimport warnings\n\nimport scipy\nfrom scipy.ndimage import zoom\nfrom skimage import measure\nfrom tensorflow.python.keras.models import Model\n\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\nimport tensorflow as tf\nimport tensorflow.keras.backend as K\n\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__)\nt_count=0\n# %% [code]\ndef 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 = '/kaggle/input/osic-pulmonary-fibrosis-progression'\nBATCH_SIZE_PRED = 1\nFEATURES = True\nADV_FEATURES = True\nC_SIGMA, C_DELTA = tf.constant(70, dtype='float32'), tf.constant(1000, dtype='float32')\nQS = [.15, .50, .85]\nIMG_SIZE = 224\nRESIZE = 224\nSEQ_LEN = 112\nQUANTIE_IMAGE_FEATURE_LEN=32\n\nLAMBDA = .8\nMDL_VERSION = 'v3'\nMODELS_PATH = '.'\n\n\nresize_dims = (112, 224, 224)\n\n\n# %% [code]\ntrain = 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')\n\n# %% [code]\ndata = 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 = [\n        'FVC',\n        'Percent',\n        'min_week_FVC'\n    ]\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)\n\n# %% [code]\nfeat_cols = [\n    'Female', 'Male',\n    'Currently smokes', '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])\nscaler = MinMaxScaler()\ndata[scale_cols] = scaler.fit_transform(data[scale_cols])\nfeat_cols.extend(scale_cols)\nprint('all data columns:', data.columns)\n\n# %% [code]\ntrain = data.loc[data.SPLIT == 'train']\ntest = data.loc[data.SPLIT == 'val']\nsubm = data.loc[data.SPLIT == 'test']\ndel data\nsubm.head()\n\nclass CropBoundingBox:\n    @staticmethod\n    def bounding_box(img3d: np.array):\n        mid_img = img3d[int(img3d.shape[0] / 2)]\n        same_first_row = (mid_img[0, :] == mid_img[0, 0]).all()\n        same_first_col = (mid_img[:, 0] == mid_img[0, 0]).all()\n        if same_first_col and same_first_row:\n            return True\n        else:\n            return False\n\n    def __call__(self, sample):\n        image, data = sample['image'], sample['metadata']\n        if not self.bounding_box(image):\n            return sample\n\n        mid_img = image[int(image.shape[0] / 2)]\n        r_min, r_max = None, None\n        c_min, c_max = None, None\n        for row in range(mid_img.shape[0]):\n            if not (mid_img[row, :] == mid_img[0, 0]).all() and r_min is None:\n                r_min = row\n            if (mid_img[row, :] == mid_img[0, 0]).all() and r_max is None \\\n                    and r_min is not None:\n                r_max = row\n                break\n\n        for col in range(mid_img.shape[1]):\n            if not (mid_img[:, col] == mid_img[0, 0]).all() and c_min is None:\n                c_min = col\n            if (mid_img[:, col] == mid_img[0, 0]).all() and c_max is None \\\n                    and c_min is not None:\n                c_max = col\n                break\n\n        image = image[:, r_min:r_max, c_min:c_max]\n        return {'image': image, 'metadata': data}\n\n\n# %% [markdown]\n# ### 3.2.2. convert_to_hu.py\n# Credits to [Guido Zuidhof's tutorial](https://www.kaggle.com/gzuidhof/full-preprocessing-tutorial).\n\n# %% [code]\nclass ConvertToHU:\n    # def __call__(self, sample):\n    #     image, data = sample['image'], sample['metadata']\n    #\n    #     img_type = data.ImageType\n    #     is_hu = img_type[0] == 'ORIGINAL' and not (img_type[2] == 'LOCALIZER')\n    #     if not is_hu:\n    #         warnings.warn(f'Patient {data.PatientID} CT Scan not cannot be'\n    #                       f'converted to Hounsfield Units (HU).')\n    #\n    #     intercept = data.RescaleIntercept\n    #     slope = data.RescaleSlope\n    #     image = (image * slope + intercept).astype(np.int16)\n    #     return {\n    #         'features': sample['features'],\n    #         'image': image,\n    #         'metadata': data,\n    #         'target': sample['target']\n    #     }\n    def __call__(self, sample):\n        image, data = sample['image'], sample['metadata']\n\n        image = image.astype(np.int16)\n        image[image == -2000] = 0\n        intercept = data.RescaleIntercept\n        slope = data.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 {\n\n            'image': image,\n            'metadata': data,\n\n        }\n\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\n# %% [markdown]\n# ### 3.2.3. resize.py\n\n# %% [code]\nclass Resize:\n    def __init__(self, output_size):\n        assert isinstance(output_size, tuple)\n        self.output_size = output_size\n\n    def __call__(self, sample):\n        image, data = sample['image'], sample['metadata']\n        resize_factor = np.array(self.output_size) / np.array(image.shape)\n        image = zoom(image, resize_factor, mode='nearest')\n        # plot_3d(image, 0)\n        return {'image': image, 'metadata': data}\n\n\nclass ReSample:\n\n    def __call__(self, sample):\n        image, data = sample['image'], sample['metadata']\n\n        image = self.resample(image, data)\n        return {'image': image, 'metadata': data}\n\n    def resample(self,image_data, scan, new_spacing=[1, 1, 1]):\n        # Determine current pixel spacing\n\n        # spacing = map(float, ([scan.SliceThickness] + scan.PixelSpacing))\n        # spacing = np.array(list(spacing))\n        spacing = np.array([scan.SliceThickness] + list(scan.PixelSpacing), dtype=np.float32)\n#         print(\"old space\\t\",spacing)\n#         print(\"new space\\t\",new_spacing)\n        resize_factor = spacing / new_spacing\n        new_real_shape = image_data.shape * resize_factor\n        new_shape = np.round(new_real_shape)\n        real_resize_factor = new_shape / image_data.shape\n        new_spacing = spacing / real_resize_factor\n#         print(\"Shape before resampling\\t\", image_data.shape)\n        image_data = scipy.ndimage.interpolation.zoom(image_data, real_resize_factor, mode='nearest')\n#         print(\"Shape after resampling\\t\", image_data.shape)\n        return image_data\n# %% [code]\nclass DataGenOsic(Sequence):\n    def __init__(self, df, image_features,tab_cols,\n                 batch_size=8, mode='fit', shuffle=False,\n                 aug=None, resize=None, 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.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        self.image_features = image_features\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), 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\n\n            imgs_seq = self.get_imgs_feature(pat_id)\n            # load the pretrained model\n\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\n\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.resize:\n    #             img = cv2.resize(img, (self.resize, self.resize))\n    #         img = img.astype(np.float32)\n    #         img = (img - (-1000)) / 1200\n    #         img = np.repeat(img[..., np.newaxis], 1, -1)\n    #         seq_imgs.append(img)\n    #     return np.array(seq_imgs).astype(np.float32)\n    def get_imgs_seq(self, pat_id):\n\n\n\n        path = '/media/feng/data/osic-pulmonary-fibrosis-progression/osic-image-np-preprocess-cached-dataset-size-112'\n        image = np.load(path + '/' + pat_id+'.npy')\n        # image = image[:32, :, :]\n        image = image.astype(np.float)\n        # Normalize the image data\n        MIN_BOUND = -1000.0\n        MAX_BOUND = 400.0\n        image = (image - MIN_BOUND) / (MAX_BOUND - MIN_BOUND)\n        image[image > 1] = 1.\n        image[image < 0] = 0.\n\n        # Zero center the data\n        PIXEL_MEAN = 0.02865046213070556\n        image = image - PIXEL_MEAN\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.resize:\n        #         img = cv2.resize(img, (self.resize, self.resize))\n        #     img = img.astype(np.float32)\n        #     img = (img - (-1000)) / 1200\n        image = np.repeat(image[..., np.newaxis], 1, -1)\n        #     seq_imgs.append(img)\n        return image\n    def get_imgs_feature(self, pat_id):\n        imgs_feature = self.image_features.values\n        feature = imgs_feature[imgs_feature[:,0]==pat_id,1:]\n                # imgs_feature.ix[[0, 2, 4, 5, 7], ['Name', 'Height', 'Weight']].head()\n\n        return feature\n\n\n    def load_scan(self,pat_id):\n        global t_count\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        try:\n            # slices = [pydicom.read_file(p) for p in path.glob('*.dcm')]\n            # print(slices[0])\n            # if slices[0].ManufacturerModelName == 'OsiriX':\n            #\n            #     image = np.stack([s.pixel_array.astype(float) for s in slices])\n            #     return image, slices[0]\n            # path = \"/media/feng/data/osic-pulmonary-fibrosis-progression/train/ID00426637202313170790466\"\n\n            file_names = sorted(os.listdir(path), key=lambda x: int(os.path.splitext(x)[0]))\n#             print(file_names)\n            slices = [pydicom.read_file(str(path) + \"/\" + file_names[i]) for i in range(len(file_names))]\n\n            # show([s.pixel_array.astype(float) for s in slices][:10])\n#             print(\"files count: {}\".format(len(slices)))\n            # print(slices[0])\n            # skip files with no SliceLocation (eg scout views)\n            # slices = []\n\n            # for f in files:\n            #     if hasattr(f,'SliceThickness'):\n            #         print(\"slices[0].SliceThickness=\".format(slices[0].SliceThickness))\n            #         slices.append(f)\n            #     elif hasattr(f, 'SliceLocation'):\n            #         slices.append(f)\n            #     else:\n            #         skipcount = skipcount + 1\n            #\n            # print(\"skipped, no SliceLocation: {}\".format(skipcount))\n\n            if hasattr(slices[0], 'ImagePositionPatient'):\n#                 print(\"change the sorting with ImagePositionPatient \")\n                slices.sort(key=lambda x: float(x.ImagePositionPatient[2]))\n                slice_thickness = np.abs(\n                    float(slices[2].ImagePositionPatient[2]) - float(slices[3].ImagePositionPatient[2]))\n                if slice_thickness == 0.0:\n                    slice_thickness = slices[0].SliceThickness\n#                 print(\"Slice Thickness: %f\" % slices[0].SliceThickness)\n                # print(\"slices[0].SpacingBetweenSlices=: %f\" % slice_thickness)\n                print(\"ImagePositionPatient slice_thickness={}\".format(slice_thickness))\n                for s in slices:\n                    s.SliceThickness = slice_thickness\n            else:\n                t_count = t_count + 1\n\n        except:\n\n            # slice_thickness = np.abs(slices[2].SliceLocation - slices[3].SliceLocation)\n            #  # print(\"SliceLocation slice_thickness={}\".format(slice_thickness))\n            #  print()\n            # # except:\n            # #     print(\"No location value\")\n            # #     print(slices[0])\n            # #     slice_thickness = slices[0].SliceThickness\n\n            warnings.warn(f'Patient {slices[0].PatientID} CT scan does not '\n                          f'in the right scan order.')\n\n        # for s in slices:\n        #     s.SliceThickness = slice_thickness\n#         print(\"No_ImagePositionPatient_count:\" + str(t_count))\n        image = np.stack([s.pixel_array.astype(float) for s in slices])\n\n        # show([image])\n        return image, slices[0]\n\nclass Image_DataGenOsic(Sequence):\n    def __init__(self, df, transform,tab_cols,\n                 batch_size=8, mode='fit', shuffle=False,\n                 aug=None, resize=None, 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.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        self.transform = transform\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, 1), dtype=np.float32)\n        pats_batch = self.df[index * self.batch_size: (index + 1) * self.batch_size]\n        for i, pat_id in enumerate(pats_batch):\n            image, metadata = self.load_scan(pat_id)\n            sample = {'image': image, 'metadata': metadata}\n            # preprocess data\n            if self.transform:\n                for t in self.transform:\n                    sample = t(sample)\n            image, data = sample['image'], sample['metadata']\n\n            image = image.astype(np.float)\n            # Normalize the image data\n            MIN_BOUND = -1000.0\n            MAX_BOUND = 400.0\n            image = (image - MIN_BOUND) / (MAX_BOUND - MIN_BOUND)\n            image[image > 1] = 1.\n            image[image < 0] = 0.\n\n            # Zero center the data\n            PIXEL_MEAN = 0.02865046213070556\n            image = image - PIXEL_MEAN\n\n            image = np.repeat(image[..., np.newaxis], 1, -1)\n\n            # imgs_seq = self.get_imgs_feature(pat_id)\n\n            X_img[i, ] = image\n\n            return (X_img, X_img),\n        else:\n            raise AttributeError('mode parameter error')\n\n\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.resize:\n    #             img = cv2.resize(img, (self.resize, self.resize))\n    #         img = img.astype(np.float32)\n    #         img = (img - (-1000)) / 1200\n    #         img = np.repeat(img[..., np.newaxis], 1, -1)\n    #         seq_imgs.append(img)\n    #     return np.array(seq_imgs).astype(np.float32)\n    def get_imgs_seq(self, pat_id):\n\n\n\n        path = '/media/feng/data/osic-pulmonary-fibrosis-progression/osic-image-np-preprocess-cached-dataset-size-112'\n        image = np.load(path + '/' + pat_id+'.npy')\n        # image = image[:32, :, :]\n        image = image.astype(np.float)\n        # Normalize the image data\n        MIN_BOUND = -1000.0\n        MAX_BOUND = 400.0\n        image = (image - MIN_BOUND) / (MAX_BOUND - MIN_BOUND)\n        image[image > 1] = 1.\n        image[image < 0] = 0.\n\n        # Zero center the data\n        PIXEL_MEAN = 0.02865046213070556\n        image = image - PIXEL_MEAN\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.resize:\n        #         img = cv2.resize(img, (self.resize, self.resize))\n        #     img = img.astype(np.float32)\n        #     img = (img - (-1000)) / 1200\n        image = np.repeat(image[..., np.newaxis], 1, -1)\n        #     seq_imgs.append(img)\n        return image\n    def get_imgs_feature(self, pat_id):\n        imgs_feature = pd.read_pickle('df_predit_image_features.pkl').values\n        feature = imgs_feature[imgs_feature[:,0]==pat_id+'.npy',1:]\n                # imgs_feature.ix[[0, 2, 4, 5, 7], ['Name', 'Height', 'Weight']].head()\n        return feature\n\n\n    def load_scan(self,pat_id):\n        global t_count\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        try:\n            # slices = [pydicom.read_file(p) for p in path.glob('*.dcm')]\n            # print(slices[0])\n            # if slices[0].ManufacturerModelName == 'OsiriX':\n            #\n            #     image = np.stack([s.pixel_array.astype(float) for s in slices])\n            #     return image, slices[0]\n            # path = \"/media/feng/data/osic-pulmonary-fibrosis-progression/train/ID00426637202313170790466\"\n\n            file_names = sorted(os.listdir(path), key=lambda x: int(os.path.splitext(x)[0]))\n#             print(file_names)\n            slices = [pydicom.read_file(str(path) + \"/\" + file_names[i]) for i in range(len(file_names))]\n\n            # show([s.pixel_array.astype(float) for s in slices][:10])\n#             print(\"files count: {}\".format(len(slices)))\n            # print(slices[0])\n            # skip files with no SliceLocation (eg scout views)\n            # slices = []\n\n            # for f in files:\n            #     if hasattr(f,'SliceThickness'):\n            #         print(\"slices[0].SliceThickness=\".format(slices[0].SliceThickness))\n            #         slices.append(f)\n            #     elif hasattr(f, 'SliceLocation'):\n            #         slices.append(f)\n            #     else:\n            #         skipcount = skipcount + 1\n            #\n            # print(\"skipped, no SliceLocation: {}\".format(skipcount))\n\n            if hasattr(slices[0], 'ImagePositionPatient'):\n#                 print(\"change the sorting with ImagePositionPatient \")\n                slices.sort(key=lambda x: float(x.ImagePositionPatient[2]))\n                slice_thickness = np.abs(\n                    float(slices[2].ImagePositionPatient[2]) - float(slices[3].ImagePositionPatient[2]))\n                if slice_thickness == 0.0:\n                    slice_thickness = slices[0].SliceThickness\n#                 print(\"Slice Thickness: %f\" % slices[0].SliceThickness)\n                # print(\"slices[0].SpacingBetweenSlices=: %f\" % slice_thickness)\n#                 print(\"ImagePositionPatient slice_thickness={}\".format(slice_thickness))\n                for s in slices:\n                    s.SliceThickness = slice_thickness\n            else:\n                t_count = t_count + 1\n\n        except:\n\n            # slice_thickness = np.abs(slices[2].SliceLocation - slices[3].SliceLocation)\n            #  # print(\"SliceLocation slice_thickness={}\".format(slice_thickness))\n            #  print()\n            # # except:\n            # #     print(\"No location value\")\n            # #     print(slices[0])\n            # #     slice_thickness = slices[0].SliceThickness\n\n            warnings.warn(f'Patient {slices[0].PatientID} CT scan does not '\n                          f'in the right scan order.')\n\n        # for s in slices:\n        #     s.SliceThickness = slice_thickness\n#         print(\"No_ImagePositionPatient_count:\" + str(t_count))\n        image = np.stack([s.pixel_array.astype(float) for s in slices])\n\n        # show([image])\n        return image, slices[0]\n\n\nclass Segenment_lung:\n    def __call__(self, sample):\n        image, data = sample['image'], sample['metadata']\n        image = self.segment_lung_mask(image, fill_lung_structures=True)\n        # plot_3d(image, 0)\n        return {'image': image, 'metadata': data}\n\n\n\n    def largest_label_volume(self,im, bg=-1):\n        vals, counts = np.unique(im, return_counts=True)\n\n        counts = counts[vals != bg]\n        vals = vals[vals != bg]\n\n        if len(counts) > 0:\n            return vals[np.argmax(counts)]\n        else:\n            return None\n\n    def segment_lung_mask(self,image, fill_lung_structures=True):\n\n        # not actually binary, but 1 and 2.\n        # 0 is treated as background, which we do not want\n        binary_image = np.array(image > -320, dtype=np.int8)+1\n        labels = measure.label(binary_image)\n\n        # Pick the pixel in the very corner to determine which label is air.\n        #   Improvement: Pick multiple background labels from around the patient\n        #   More resistant to \"trays\" on which the patient lays cutting the air\n        #   around the person in half\n        background_label = labels[0,0,0]\n\n        #Fill the air around the person\n        binary_image[background_label == labels] = 2\n\n\n        # Method of filling the lung structures (that is superior to something like\n        # morphological closing)\n        if fill_lung_structures:\n            # For every slice we determine the largest solid structure\n            for i, axial_slice in enumerate(binary_image):\n                axial_slice = axial_slice - 1\n                labeling = measure.label(axial_slice)\n                l_max = self.largest_label_volume(labeling, bg=0)\n\n                if l_max is not None: #This slice contains some lung\n                    binary_image[i][labeling != l_max] = 1\n\n\n        binary_image -= 1 #Make the image actual binary\n        binary_image = 1-binary_image # Invert it, lungs are now 1\n\n        # Remove other air pockets insided body\n        labels = measure.label(binary_image, background=0)\n        l_max = self.largest_label_volume(labels, bg=0)\n        if l_max is not None: # There are air pockets\n            binary_image[labels != l_max] = 0\n\n        return binary_image\n# %% [code]\ndef 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]\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) * sq2 + 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\n\n# %% [code]\n# %%time\nmodel_file = '/kaggle/input/model-file/model_v3.h5'\nimage_model_file='/kaggle/input/image-model/model_3d_v0.h5'\n\n\ntest_image_path = f'{DATA_PATH}/test'\ntest_data = sorted(os.listdir(test_image_path), key=lambda x: os.path.splitext(x)[0])\n\ntest_dataset = Image_DataGenOsic(\n    df=test_data,\n    transform=[\n        CropBoundingBox(),\n        ConvertToHU(),\n        ReSample(),\n        # Clip(bounds=clip_bounds),\n        Resize(resize_dims),\n        Segenment_lung()],\n\n    tab_cols=feat_cols,\n    batch_size=BATCH_SIZE_PRED,\n    mode='predict',\n    shuffle=False,\n    aug=None,\n    resize=RESIZE,\n    seq_len=SEQ_LEN,\n    img_size=IMG_SIZE\n\n\n)\n\n# load the pretrained model\nimage_model = load_model(image_model_file)\nmean_model = Model(inputs=image_model.inputs, outputs=image_model.get_layer('Dec_VAE_VDraw_Mean').output)\n\npredit_image_features =mean_model.predict(test_dataset,verbose=1)\ntest_data = pd.DataFrame(test_data)\ndf_predit_image_features = pd.DataFrame(predit_image_features)\n\ndf_predit_image_features = pd.concat([test_data,df_predit_image_features],axis=1)\n\nprint(df_predit_image_features)\n# %% [code]\nsubm_datagen = DataGenOsic(\n    df=subm,\n    image_features=df_predit_image_features,\n    tab_cols=feat_cols,\n    batch_size=BATCH_SIZE_PRED,\n    mode='predict',\n    shuffle=False,\n    aug=None,\n    resize=RESIZE,\n    seq_len=QUANTIE_IMAGE_FEATURE_LEN,\n    img_size=IMG_SIZE\n)\n\nmodel = load_model(model_file, custom_objects={'qloss': qloss, 'loss': mloss(LAMBDA), 'score': score})\nprint('model loaded:', model_file)\nprint(model.summary())\n\n# %% [code]\n(Xt_img, Xt_tab), _ = subm_datagen.__getitem__(0)\nprint('val X img: ', Xt_img.shape)\nprint('val X tab: ', Xt_tab.shape)\n# fig, axes = plt.subplots(figsize=(10, 4), nrows=BATCH_SIZE_PRED, ncols=SEQ_LEN)\n# for 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}')\n# plt.show()\n\n# %% [code]\n# %%time\npreds_subm = model.predict(subm_datagen, verbose=1)\nprint('predictions shape:', preds_subm.shape)\nprint('predictions sample:', preds_subm[0])\n\n# # %% [code]\n# subm['FVC'] = preds_subm[:, 1]\n# subm['Confidence'] = preds_subm[:, 2] - preds_subm[:, 0]\n# subm[['Patient_Week','FVC','Confidence']].to_csv('submission.csv', index=False)\n# subm[['Patient_Week','FVC','Confidence']].describe().T\n\n# %% [code]\nsubm['FVC1'] = preds_subm[:, 1]\nsubm['Confidence1'] = preds_subm[:,2] - preds_subm[:,0]\n\n# get rid of unused data and show some non-empty data\nsubmission = subm[['Patient_Week','FVC','Confidence','FVC1','Confidence1']].copy()\nsubmission.loc[~submission.FVC1.isnull()].head(10)\n\n# %% [code]\nsubmission.loc[~submission.FVC1.isnull(),'FVC'] = submission.loc[~submission.FVC1.isnull(),'FVC1']\n\n\nsubmission.loc[~submission.FVC1.isnull(),'Confidence'] = submission.loc[~submission.FVC1.isnull(),'Confidence1']\n\norg_test = pd.read_csv('/kaggle/input/osic-pulmonary-fibrosis-progression/test.csv')\n\nfor i in range(len(org_test)):\n    submission.loc[submission['Patient_Week']==org_test.Patient[i]+'_'+str(org_test.Weeks[i]), 'FVC'] = org_test.FVC[i]\n    submission.loc[submission['Patient_Week']==org_test.Patient[i]+'_'+str(org_test.Weeks[i]), 'Confidence'] = 70\n\nsubmission[[\"Patient_Week\",\"FVC\",\"Confidence\"]].to_csv(\"submission.csv\", index = False)","execution_count":null,"outputs":[]},{"metadata":{"_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","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}