{"cells":[{"metadata":{},"cell_type":"markdown","source":"This notebook was built by assembling pieces of other kaggle notebook in this competition"},{"metadata":{},"cell_type":"markdown","source":"# About Pulmonary Fibrosis\nreference: https://www.kaggle.com/awwaldiekaramapepple/your-starter-notebook-for-osic/edit"},{"metadata":{},"cell_type":"markdown","source":"![](https://www.osicild.org/uploads/1/2/2/7/122798879/editor/kaggle-v01-clipped_2.png?1569348761)"},{"metadata":{},"cell_type":"markdown","source":"> **What do you need to do?**\n\n> In this competition, you’ll predict a patient’s severity of decline in lung function based on a CT scan of their lungs. You’ll determine lung function based on output from a spirometer, which measures the volume of air inhaled and exhaled. The challenge is to use machine learning techniques to make a prediction with the image, metadata, and baseline FVC as input.\n\n\n> **Submission**\n\n> The test set consists of Three_Patient_Week(s) per patient. You need to predict the **Forced vital capacity(FVC)** i.e. volume of air exhaled and the **confidence** value in your prediction.\n\n\n>**Your Evaluation Metric**\n\n> For each true FVC measurement, you will predict both an FVC and a confidence measure (standard deviation σ). The metric is computed as:\n\n\\begin{equation} \n\\sigma_{clipped} = max(\\sigma, 70)\\\\\n\\Delta = min ( |FVC_{true} - FVC_{predicted}|, 1000 )\\\\\nmetric = -   \\frac{\\sqrt{2} \\Delta}{\\sigma_{clipped}} - \\ln ( \\sqrt{2} \\sigma_{clipped} )\n\\end{equation}\n\n\n> A great notebook to understand the metric --> https://www.kaggle.com/rohanrao/osic-understanding-laplace-log-likelihood"},{"metadata":{"trusted":true},"cell_type":"code","source":"import pandas as pd\nimport numpy as np\nimport matplotlib.pyplot as plt\nimport pydicom\n\nimport os\nimport time","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"import imageio\nfrom IPython.display import Image\nfrom skimage.measure import label, regionprops","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"INPUT_FOLDER = '/kaggle/input/osic-pulmonary-fibrosis-progression/train/'\n\npatients = os.listdir(INPUT_FOLDER)\npatients.sort()\n\nprint('some patients: \\n', '\\n'.join(patients[:5]))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def load_scan(path):\n    '''\n    Loads scans from folder into a list\n    Args path of images to load\n    Returns images path in a list\n    '''\n    slices = [pydicom.read_file(path +'/'+ s) for s in os.listdir(path)]\n    slices.sort(key = lambda x: int(x.InstanceNumber))\n    \n    try:\n        slice_thickness = np.abs(slices[0].ImagePositionPatient[2] - slices[0].ImagePositionPatient[2])\n    except:\n        slice[0].slice_location - slice[1].slice_location\n    for s in slices:\n        s.SliceThickness = slice_thickness\n    \n    return slices\n        ","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"the unit of measurement of CT scans is hounsfield unit (***hu***)\nformulae for hounsfield unit\n* hu = m* P + b\n* m = rescale slope\n* p = rescale intercept\n* b = pixel array"},{"metadata":{},"cell_type":"markdown","source":"REFEERENCE KERNELS\n%\n* https://www.kaggle.com/aadhavvignesh/lung-segmentation-by-marker-controlled-watershed#Conclusion:\n\n* https://www.raddq.com/dicom-processing-segmentation-visualization-in-python"},{"metadata":{"trusted":true},"cell_type":"code","source":"#helper function to convert dicom images to houndsfield\n\ndef get_pixels_hu(scans):\n    '''\n    converts raw image files into hounsfield unit\n    Arguments: raw images\n    Returns: images numpy array\n    '''\n    image = np.stack([s.pixel_array for s in scans])\n    image = image.astype(np.int16)\n    # Since the scanning equipment is cylindrical in nature and image output is square,\n    # we set the out-of-scan pixels to 0\n    \n    image[image ==-2000] = 0\n    intercept = scans[0].RescaleIntercept\n    slope = scans[0].RescaleSlope\n    \n    if slope != 1:\n        image = slope * image.astype(np.float64)\n        image = image.astype(np.int16)\n        \n    image += np.int16(intercept)\n    \n    return np.array(image, dtype=np.int16)\n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"path = INPUT_FOLDER + patients[24]\ntest_patient_scan = load_scan(path)\ntest_patient_images = get_pixels_hu(test_patient_scan)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"path = INPUT_FOLDER + patients[24]\nslices = [pydicom.read_file(path +'/'+ s) for s in os.listdir(path)]","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"we'll take random slices and perform segmentation"},{"metadata":{"trusted":true},"cell_type":"code","source":"plt.imshow(test_patient_images[12]) ;\nplt.title('original image slice 12');","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"creating markers"},{"metadata":{"trusted":true},"cell_type":"code","source":"from skimage import measure, morphology, segmentation\nimport scipy.ndimage as ndimage","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def generate_markers(image):\n    '''\n    Generate markers for a given image\n    Arguments: image\n    returns : Internal marker, external marker, watershed marker\n    '''\n    #creation the internal marker\n    marker_internal = image < -400\n    marker_internal = segmentation.clear_border(marker_internal)\n    marker_internal_labels = measure.label(marker_internal)\n    \n    areas = [r.area for r in measure.regionprops(marker_internal_labels)]\n    areas.sort()\n    \n    if len(areas) > 2:\n        for region in measure.regionprops(marker_internal_labels):\n            if region.area < areas[-2]:\n                for coordinates in region.coords:\n                    marker_internal_labels[coordinates[0], coordinates[1]] == 0\n                    \n    marker_internal = marker_internal_labels > 0\n    \n    #creation of external marker\n    external_a = ndimage.binary_dilation(marker_internal, iterations=10)\n    external_b = ndimage.binary_dilation(marker_internal, iterations=55)\n    marker_external = external_b ^ external_a\n    \n    #creation of watershed marker\n    marker_watershed = np.zeros((512,512), dtype=np.int)\n    marker_watershed += marker_internal * 255\n    marker_watershed += marker_external * 128\n    \n    return marker_internal, marker_external, marker_watershed\n    ","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"test_patient_internal, test_patient_external, test_patient_watershed = generate_markers(\n                                                                        test_patient_images[12])\nf, (ax1, ax2, ax3) = plt.subplots(1, 3, sharey=True, figsize= (15,15))\n\nax1.imshow(test_patient_internal)#, cmap='gray')\nax1.set_title('internal_marker')\nax1.axis('off')\n\nax2.imshow(test_patient_external)\nax2.set_title('external_marker')\n\nax3.imshow(test_patient_watershed)\nax3.set_title('watershed')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"f, (ax1, ax2, ax3) = plt.subplots(1, 3, sharey=True, figsize=(15,15))\n\nax1.imshow(test_patient_internal, cmap='gray')\nax1.set_title(\"Internal Marker\")\nax1.axis('off')\n\nax2.imshow(test_patient_external, cmap='gray')\nax2.set_title(\"External Marker\")\nax2.axis('off')\n\nax3.imshow(test_patient_watershed, cmap='gray')\nax3.set_title(\"Watershed Marker\")\nax3.axis('off')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"#list to store computation times and iterations\ncompute_time = []\niter_titles = []","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"The most basic morphological operations are dilation and erosion. Dilation adds pixels to the boundaries of objects in an image, while erosion removes pixels on object boundaries. The number of pixels added or removed from the objects in an image depends on the size and shape of the structuring element used to process the image. In the morphological dilation and erosion operations, the state of any given pixel in the output image is determined by applying a rule to the corresponding pixel and its neighbors in the input image. The rule used to process the pixels defines the operation as a dilation or an erosion."},{"metadata":{"trusted":true},"cell_type":"code","source":"def seperate_lungs(image, iterations=1):\n    \"\"\"\n    Segments lungs using various techniques\n    Parameters: image (scan images) iteration (number of iteration)\n    returns:  -Segmented Lung, -Lung Filter, -Outline Lung, \n              -watershed Lung, -Sobel Gradient\n    \"\"\"\n    # Store the start time\n    start = time.time()\n    \n    marker_internal, marker_external, marker_watershed = generate_markers(image)\n    \n    \"\"\"\n    creation of sobel gradient\n    \"\"\"\n    # Sobel Gradient\n    sobel_filtered_dx = ndimage.sobel(image, 1)\n    sobel_filtered_dy = ndimage.sobel(image, 0)\n    sobel_gradient = np.hypot(sobel_filtered_dx, sobel_filtered_dy)\n    sobel_gradient *= 255.0 / np.max(sobel_gradient)\n    \n    \"\"\"\n    using algorithm watershed\n    \n    we pass the image convoluted by sobel and watershed marker\n    to morphology watershed and get a matrix matrix labelled using \n    the watershed segmentation algorithm\n    \"\"\"\n    watershed = morphology.watershed(sobel_gradient, marker_watershed)\n    \"\"\"\n    Reducing the image to outlines after watershed algorithm\n    \"\"\"\n    outline = ndimage.morphological_gradient(watershed, size=(3,3))\n    outline = outline.astype(bool)\n    \"\"\"\n    Black Top-Hart morphology:\n    \n    \n    \"\"\"\n    #structuring element used for filter\n    blackhat_struct = [[0, 0, 1, 1, 1, 0, 0],\n                        [0, 1, 1, 1, 1, 1, 0],\n                        [1, 1, 1, 1, 1, 1, 1],\n                        [1, 1, 1, 1, 1, 1, 1],\n                        [1, 1, 1, 1, 1, 1, 1],\n                        [0, 1, 1, 1, 1, 1, 0],\n                        [0, 0, 1, 1, 1, 0, 0]]\n    blackhat_struct = ndimage.iterate_structure(blackhat_struct, iterations)\n    #Perform black top-hat filter\n    outline = ndimage.black_tophat(outline, structure= blackhat_struct)\n    \"\"\"\n    Generate internal filter using internal marker and outline\n    \"\"\"\n    lung_filter = np.bitwise_or(marker_internal, outline)\n    lung_filter = ndimage.morphology.binary_closing(lung_filter, structure=np.ones((3,3)), iterations=1)\n    \"\"\"\n    Segment lung using lungfilter and the image\n    \"\"\"\n    segmented = np.where(lung_filter==1, image, -2000* np.ones((512, 512)))\n    \n    #Append Computation time\n    end = time.time()\n    compute_time.append(end - start)\n    iter_titles.append(\"{num} iterations\".format(num=iterations))\n    \n    return segmented, lung_filter, outline, watershed, sobel_gradient","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"Comparison of iterations with time\nWe'll be checking for iterations in the range of 1-8. iterations = 1 is the default for the seperate_lungs function."},{"metadata":{"trusted":true},"cell_type":"code","source":"for itr in range(1,9):\n    (test_segmented, test_lung_filter, test_outline,\n    test_watershed, test_sobel_gradient) = seperate_lungs(test_patient_images[12], itr)#test_patient_images[12]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"itr_dict = {'Iterations' : iter_titles, 'computation time (in seconds)' : compute_time}\ncolors = ['#30336b',] * 8\ncolors[0] = '#ed4d4b'\n\nimport plotly.express as px\nimport plotly.graph_objects as go\n\nfig = go.Figure(data=[go.Bar(\n            x= itr_dict['Iterations'],\n            y= itr_dict['computation time (in seconds)'],\n            marker_color = colors\n            )])\n\nfig.update_traces(texttemplate= '%{y:.3s}', textposition= 'outside')\n\nfig.update_layout(\n            title = 'Iterations vs computation times',\n            yaxis = dict(\n                    title='Computation time (in seconds)',\n                    titlefont_size= 16,\n                    tickfont_size= 14,\n                        ),\n            autosize= False,\n            width= 800,\n            height= 700,)\n\nfig.show()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"f, ax = plt.subplots(1,2, sharey=True, figsize=(12,12))\nax[0].imshow(test_sobel_gradient)\nax[0].set_title('test sobel gradient')\nax[0].axis('off')\n\nax[1].imshow(test_watershed)\nax[1].set_title('test_watershed')\nax[1].axis('off')\n\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"f, (ax1, ax2, ax3) = plt.subplots(1,3, sharex=True, figsize=(12,12))\n\nax1.imshow(test_outline)\nax1.set_title('test_outline')\nax2.imshow(test_lung_filter)\nax2.set_title('test_lung_filter')\nax2.axis('off')\nax3.imshow(test_segmented)\nax3.set_title('test_segment')\nax3.axis('off')\n\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"f, ax = plt.subplots(1,2, sharey=True, figsize=(14, 12))\nax[0].imshow(test_patient_images[12])\nax[0].set_title('original image')\nax[0].axis('off')\n\nax[1].imshow(test_segmented)\nax[1].set_title('segmented image')\nax[1].axis('off')\n\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def set_lungwin(img, hu=[-1200., 600.]):\n    lungwin = np.array(hu)\n    newimg = (img-lungwin[0]) / (lungwin[1]-lungwin[0])\n    newimg[newimg < 0] = 0\n    newimg[newimg > 1] = 1\n    newimg = (newimg * 255).astype('uint8')\n    return newimg\n\n\nscans = load_scan('../input/osic-pulmonary-fibrosis-progression/train/ID00007637202177411956430/')\nscan_array = set_lungwin(get_pixels_hu(scans))\n\nimageio.mimsave(\"/tmp/gif.gif\", scan_array, duration=0.0001)\nImage(filename=\"/tmp/gif.gif\", format='png')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"slices = [path +'/'+ s for s in os.listdir(path)]\nslices[5]","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"#sample_image = pydicom.dcmread(scans[7])\nsample_image = pydicom.dcmread(slices[7])\nimg = sample_image.pixel_array\n\nplt.imshow(img, cmap='gray') ;\n#print(img.value)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"img = (img + sample_image.RescaleIntercept) / sample_image.RescaleSlope\nimg = img < -400 #HU unit range for lung CT scans\nf, (ax1, ax2) = plt.subplots(1,2, figsize=(10, 10))\n\nax1.imshow(img, cmap='gray')\nax1.set_title('Binary mask image')\n\nimg = segmentation.clear_border(img)\nax2.imshow(img, cmap = 'gray')\nax2.set_title('cleaned border image')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"#### Image Labelling"},{"metadata":{"trusted":true},"cell_type":"code","source":"img = label(img)\nplt.imshow(img, cmap='gray') ;","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"len([r for r in regionprops(img)])\nareas = [r.area for r in regionprops(img)]\nareas.sort()\nif len(areas) > 2:\n    for region in regionprops(img):\n        if region.area < areas[-2]:\n            for coordinates in region.coords:\n                img[coordinates[0], coordinates[1]] = 0\nimg = img > 0\nplt.imshow (img)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"https://www.kaggle.com/khoongweihao/efficientnets-quantile-regression-inference"},{"metadata":{"trusted":true},"cell_type":"code","source":"import tensorflow as tf\nimport random\n\nfrom tensorflow.keras import Model\nimport tensorflow.keras.backend as K\nimport tensorflow.keras.layers as L\nimport tensorflow.keras.models as M\nfrom tensorflow.keras.optimizers import Nadam\nfrom tensorflow.keras.utils import Sequence\n\nimport seaborn as sns\nfrom PIL import Image","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"from tqdm.notebook import tqdm\nimport cv2","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def seed_everything(seed= 2020):\n    random.seed(seed)\n    os.environ['PYTHONHASHSEED'] = str(seed)\n    np.random.seed(seed)\n    tf.random.set_seed(seed)\n    tf.random.set_seed(seed)\n    \nseed_everything(42)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"#tf.compact allows us to write code that works both in tf 1.x and 2.x e.g tf.compact.v2 allows us\n#to use things introduced in 2.x from 1.x\nconfig = tf.compat.v1.ConfigProto()\nconfig.gpu_options.allow_growth = True\nsession = tf.compat.v1.Session(config= config)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train = pd.read_csv('../input/osic-pulmonary-fibrosis-progression/train.csv')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"train.sample(5)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def get_tab(df):\n    vector = [(df.Age.values[0] - 30) / 30]\n    \n    if df.Sex.values[0] == 'male':\n        vector.append(0)\n    else:\n        vector.append(1)\n        \n    if df['SmokingStatus'].values[0] == 'Never Smoked':\n        vector.extend([0, 0])\n    elif df['SmokingStatus'].values[0] == 'Ex-smoker':\n        vector.extend([1,1])\n    elif df['SmokingStatus'].values[0] == 'Currently Smokes':\n        vector.extend([0, 1])\n    else:\n        vector.extend([1, 0])\n    \n    return np.array(vector)\n","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"\nIn summary efficientnets is a pretrained model\n* https://medium.com/analytics-vidhya/image-classification-with-efficientnet-better-performance-with-computational-efficiency-f480fdb00ac6\n\nModel scaling is about scaling the existing model in terms of model depth, model width, and less popular input image resolution to improve the performance of the model. Depth wise scaling is most popular amongst all, e.g. ResNet can be scaled from Resnet18 to ResNet200\nResNet200 delivers better performance than ResNet18, there is a problem with traditional manual scaling method, after a certain level, scaling doesn’t improve performance. It starts to affect adversely by degrading performance.\nCompound scaling instead of scaling only one model attribute out of depth, width, and resolution; strategically scaling all three of them together delivers better results\nAuthors observed that mobile scaling can be used on any CNN architecture and it works fine but the overall performance very much depends on baseline architecture. With that, they came up with the brand new base architecture and named it EfficientNet-B0."},{"metadata":{},"cell_type":"markdown","source":"## Efficientnet\n\nRefrence kernel: https://www.kaggle.com/khoongweihao/efficientnets-quantile-regression-inference/data"},{"metadata":{"trusted":true},"cell_type":"code","source":"!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":{"trusted":true},"cell_type":"code","source":"A = {} \nTAB = {} \nP = [] \nfor i, p in tqdm(enumerate(train.Patient.unique())):\n    sub = train.loc[train.Patient == p, :] \n    fvc = sub.FVC.values\n    weeks = sub.Weeks.values\n    c = np.vstack([weeks, np.ones(len(weeks))]).T\n    a, b = np.linalg.lstsq(c, fvc)[0]\n    \n    A[p] = a\n    TAB[p] = get_tab(sub)\n    P.append(p)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def get_img(path):\n    d = pydicom.dcmread(path)\n    return cv2.resize(d.pixel_array / 2**11, (512, 512))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"from tensorflow.keras.utils import Sequence\n\nclass IGenerator(Sequence):\n    BAD_ID = ['ID00011637202177653955184', 'ID00052637202186188008618']\n    def __init__(self, keys, a, tab, batch_size=32):\n        self.keys = [k for k in keys if k not in self.BAD_ID]\n        self.a = a\n        self.tab = tab\n        self.batch_size = batch_size\n        \n        self.train_data = {}\n        for p in train.Patient.values:\n            self.train_data[p] = os.listdir(f'../input/osic-pulmonary-fibrosis-progression/train/{p}/')\n    \n    def __len__(self):\n        return 1000\n    \n    def __getitem__(self, idx):\n        x = []\n        a, tab = [], [] \n        keys = np.random.choice(self.keys, size = self.batch_size)\n        for k in keys:\n            try:\n                i = np.random.choice(self.train_data[k], size=1)[0]\n                img = get_img(f'../input/osic-pulmonary-fibrosis-progression/train/{k}/{i}')\n                x.append(img)\n                a.append(self.a[k])\n                tab.append(self.tab[k])\n            except:\n                print(k, i)\n       \n        x,a,tab = np.array(x), np.array(a), np.array(tab)\n        x = np.expand_dims(x, axis=-1)\n        return [x, tab] , a","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"from tensorflow.keras.layers import (\n    Dense, Dropout, Activation, Flatten, Input, BatchNormalization, GlobalAveragePooling2D, Add, Conv2D, AveragePooling2D, \n    LeakyReLU, Concatenate \n)\nimport efficientnet.tfkeras as efn\n\ndef get_efficientnet(model, shape):\n    models_dict = {\n        'b0': efn.EfficientNetB0(input_shape=shape,weights=None,include_top=False),\n        'b1': efn.EfficientNetB1(input_shape=shape,weights=None,include_top=False),\n        'b2': efn.EfficientNetB2(input_shape=shape,weights=None,include_top=False),\n        'b3': efn.EfficientNetB3(input_shape=shape,weights=None,include_top=False),\n        'b4': efn.EfficientNetB4(input_shape=shape,weights=None,include_top=False),\n        'b5': efn.EfficientNetB5(input_shape=shape,weights=None,include_top=False),\n        'b6': efn.EfficientNetB6(input_shape=shape,weights=None,include_top=False),\n        'b7': efn.EfficientNetB7(input_shape=shape,weights=None,include_top=False)\n    }\n    return models_dict[model]\n\ndef build_model(shape=(512, 512, 1), model_class=None):\n    inp = Input(shape=shape)\n    base = get_efficientnet(model_class, shape)\n    x = base(inp)\n    x = GlobalAveragePooling2D()(x)\n    inp2 = Input(shape=(4,))\n    x2 = tf.keras.layers.GaussianNoise(0.2)(inp2)\n    x = Concatenate()([x, x2]) \n    x = Dropout(0.5)(x) \n    x = Dense(1)(x)\n    model = Model([inp, inp2] , x)\n    \n    weights = [w for w in os.listdir('../input/osic-model-weights') if model_class in w][0]\n    model.load_weights('../input/osic-model-weights/' + weights)\n    return model\n\nmodel_classes = ['b5'] #['b0','b1','b2','b3',b4','b5','b6','b7']\nmodels = [build_model(shape=(512, 512, 1), model_class=m) for m in model_classes]\nprint('Number of models: ' + str(len(models)))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"from sklearn.model_selection import train_test_split \n\ntr_p, vl_p = train_test_split(P, \n                              shuffle=True, \n                              train_size= 0.8) \n","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sns.distplot(list(A.values()))","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"def score(fvc_true, fvc_pred, sigma):\n    sigma_clip = np.maximum(sigma, 70) # changed from 70, trie 66.7 too\n    delta = np.abs(fvc_true - fvc_pred)\n    delta = np.minimum(delta, 1000)\n    sq2 = np.sqrt(2)\n    metric = (delta / sigma_clip)*sq2 + np.log(sigma_clip* sq2)\n    return np.mean(metric)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"subs = []\nfor model in models:\n    metric = []\n    for q in tqdm(range(1, 10)):\n        m = []\n        for p in vl_p:\n            x = [] \n            tab = [] \n\n            if p in ['ID00011637202177653955184', 'ID00052637202186188008618']:\n                continue\n\n            ldir = os.listdir(f'../input/osic-pulmonary-fibrosis-progression/train/{p}/')\n            for i in ldir:\n                if int(i[:-4]) / len(ldir) < 0.8 and int(i[:-4]) / len(ldir) > 0.15:\n                    x.append(get_img(f'../input/osic-pulmonary-fibrosis-progression/train/{p}/{i}')) \n                    tab.append(get_tab(train.loc[train.Patient == p, :])) \n            if len(x) < 1:\n                continue\n            tab = np.array(tab) \n\n            x = np.expand_dims(x, axis=-1) \n            _a = model.predict([x, tab]) \n            a = np.quantile(_a, q / 10)\n\n            percent_true = train.Percent.values[train.Patient == p]\n            fvc_true = train.FVC.values[train.Patient == p]\n            weeks_true = train.Weeks.values[train.Patient == p]\n\n            fvc = a * (weeks_true - weeks_true[0]) + fvc_true[0]\n            percent = percent_true[0] - a * abs(weeks_true - weeks_true[0])\n            m.append(score(fvc_true, fvc, percent))\n        print(np.mean(m))\n        metric.append(np.mean(m))\n\n    q = (np.argmin(metric) + 1)/ 10\n\n    sub = pd.read_csv('../input/osic-pulmonary-fibrosis-progression/sample_submission.csv') \n    test = pd.read_csv('../input/osic-pulmonary-fibrosis-progression/test.csv') \n    A_test, B_test, P_test,W, FVC= {}, {}, {},{},{} \n    STD, WEEK = {}, {} \n    for p in test.Patient.unique():\n        x = [] \n        tab = [] \n        ldir = os.listdir(f'../input/osic-pulmonary-fibrosis-progression/test/{p}/')\n        for i in ldir:\n            if int(i[:-4]) / len(ldir) < 0.8 and int(i[:-4]) / len(ldir) > 0.15:\n                x.append(get_img(f'../input/osic-pulmonary-fibrosis-progression/test/{p}/{i}')) \n                tab.append(get_tab(test.loc[test.Patient == p, :])) \n        if len(x) <= 1:\n            continue\n        tab = np.array(tab) \n\n        x = np.expand_dims(x, axis=-1) \n        _a = model.predict([x, tab]) \n        a = np.quantile(_a, q)\n        A_test[p] = a\n        B_test[p] = test.FVC.values[test.Patient == p] - a*test.Weeks.values[test.Patient == p]\n        P_test[p] = test.Percent.values[test.Patient == p] \n        WEEK[p] = test.Weeks.values[test.Patient == p]\n\n    for k in sub.Patient_Week.values:\n        p, w = k.split('_')\n        w = int(w) \n\n        fvc = A_test[p] * w + B_test[p]\n        sub.loc[sub.Patient_Week == k, 'FVC'] = fvc\n        sub.loc[sub.Patient_Week == k, 'Confidence'] = (\n            P_test[p] - A_test[p] * abs(WEEK[p] - w) \n    ) \n\n    _sub = sub[[\"Patient_Week\",\"FVC\",\"Confidence\"]].copy()\n    subs.append(_sub)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"N = len(subs)\nsub = subs[0].copy() # ref\nsub[\"FVC\"] = 0\nsub[\"Confidence\"] = 0\nfor i in range(N):\n    sub[\"FVC\"] += subs[0][\"FVC\"] * (1/N)\n    sub[\"Confidence\"] += subs[0][\"Confidence\"] * (1/N)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sub.head()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sub.to_csv('submission.csv', index= False)","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}