{"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\n\nimport pydicom\nimport pandas as pd\nimport numpy as np \nimport tensorflow as tf \nimport matplotlib.pyplot as plt \n\nfrom tqdm.notebook import tqdm ","metadata":{"_cell_guid":"79c7e3d0-c299-4dcb-8224-4455121ee9b0","_uuid":"d629ff2d2480ee46fbb7e2d37f6b5fab8052498a","execution":{"iopub.status.busy":"2021-11-27T06:27:42.279653Z","iopub.execute_input":"2021-11-27T06:27:42.280027Z","iopub.status.idle":"2021-11-27T06:27:47.141725Z","shell.execute_reply.started":"2021-11-27T06:27:42.279991Z","shell.execute_reply":"2021-11-27T06:27:47.140912Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Decay theory\nInput for test:\n   * FVC in n week\n   * Percent in n week \n   * Age\n   * Sex\n   * Smoking status\n   * CT in n week\n   \nResult:\n   * FVC in any week\n   * percent in any week\n   \n$FVC = a.quantile(0.75) * (week - week_{test}) + FVC_{test}$\n\n$Confidence = Percent + a.quantile(0.75) * abs(week - week_{test}) $\n\nSo let's try predict coefficient a. ","metadata":{}},{"cell_type":"code","source":"train = pd.read_csv('../input/osic-pulmonary-fibrosis-progression/train.csv') ","metadata":{"execution":{"iopub.status.busy":"2021-11-27T06:27:47.14346Z","iopub.execute_input":"2021-11-27T06:27:47.143775Z","iopub.status.idle":"2021-11-27T06:27:47.16348Z","shell.execute_reply.started":"2021-11-27T06:27:47.143748Z","shell.execute_reply":"2021-11-27T06:27:47.162844Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train.head()","metadata":{"execution":{"iopub.status.busy":"2021-11-27T06:27:47.164924Z","iopub.execute_input":"2021-11-27T06:27:47.165417Z","iopub.status.idle":"2021-11-27T06:27:47.18649Z","shell.execute_reply.started":"2021-11-27T06:27:47.16538Z","shell.execute_reply":"2021-11-27T06:27:47.185718Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train.SmokingStatus.unique()","metadata":{"execution":{"iopub.status.busy":"2021-11-27T06:27:47.187879Z","iopub.execute_input":"2021-11-27T06:27:47.188267Z","iopub.status.idle":"2021-11-27T06:27:47.200108Z","shell.execute_reply.started":"2021-11-27T06:27:47.18823Z","shell.execute_reply":"2021-11-27T06:27:47.199385Z"},"trusted":true},"execution_count":null,"outputs":[]},{"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    return np.array(vector) ","metadata":{"execution":{"iopub.status.busy":"2021-11-27T06:27:47.204703Z","iopub.execute_input":"2021-11-27T06:27:47.204978Z","iopub.status.idle":"2021-11-27T06:27:47.215098Z","shell.execute_reply.started":"2021-11-27T06:27:47.204951Z","shell.execute_reply":"2021-11-27T06:27:47.213915Z"},"trusted":true},"execution_count":null,"outputs":[]},{"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)","metadata":{"execution":{"iopub.status.busy":"2021-11-27T06:27:47.218781Z","iopub.execute_input":"2021-11-27T06:27:47.219477Z","iopub.status.idle":"2021-11-27T06:27:47.563705Z","shell.execute_reply.started":"2021-11-27T06:27:47.219437Z","shell.execute_reply":"2021-11-27T06:27:47.562642Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## CNN for coeff prediction","metadata":{}},{"cell_type":"code","source":"def get_img(path):\n    d = pydicom.dcmread(path)\n    return cv2.resize(d.pixel_array / 2**11, (512, 512))","metadata":{"execution":{"iopub.status.busy":"2021-11-27T06:27:47.565081Z","iopub.execute_input":"2021-11-27T06:27:47.565439Z","iopub.status.idle":"2021-11-27T06:27:47.571981Z","shell.execute_reply.started":"2021-11-27T06:27:47.5654Z","shell.execute_reply":"2021-11-27T06:27:47.570875Z"},"trusted":true},"execution_count":null,"outputs":[]},{"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","metadata":{"execution":{"iopub.status.busy":"2021-11-27T06:27:47.573794Z","iopub.execute_input":"2021-11-27T06:27:47.574439Z","iopub.status.idle":"2021-11-27T06:27:47.591048Z","shell.execute_reply.started":"2021-11-27T06:27:47.574313Z","shell.execute_reply":"2021-11-27T06:27:47.589902Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from tensorflow.keras.layers import (\n    Dense, Dropout, Activation, Flatten, Input, BatchNormalization, GlobalAveragePooling2D, Add, Conv2D, AveragePooling2D, \n    LeakyReLU, Concatenate \n)\n\nfrom tensorflow.keras import Model\nfrom tensorflow.keras.optimizers import Nadam\n\ndef get_model(shape=(512, 512, 1)):\n    def res_block(x, n_features):\n        _x = x\n        x = BatchNormalization()(x)\n        x = LeakyReLU(0.05)(x)\n    \n        x = Conv2D(n_features, kernel_size=(3, 3), strides=(1, 1), padding='same')(x)\n        x = Add()([_x, x])\n        return x\n    \n    inp = Input(shape=shape)\n    \n    # 512\n    x = Conv2D(32, kernel_size=(3, 3), strides=(1, 1), padding='same')(inp)\n    x = BatchNormalization()(x)\n    x = LeakyReLU(0.05)(x)\n    \n    x = Conv2D(32, kernel_size=(3, 3), strides=(1, 1), padding='same')(x)\n    x = BatchNormalization()(x)\n    x = LeakyReLU(0.05)(x)\n    \n    x = AveragePooling2D(pool_size=(2, 2), strides=(2, 2))(x)\n    \n    # 256\n    x = Conv2D(8, kernel_size=(3, 3), strides=(1, 1), padding='same')(x)\n    for _ in range(2):\n        x = res_block(x, 8)\n    x = AveragePooling2D(pool_size=(2, 2), strides=(2, 2))(x)\n    \n    # 128\n    x = Conv2D(16, kernel_size=(3, 3), strides=(1, 1), padding='same')(x)\n    for _ in range(2):\n        x = res_block(x, 16)\n    x = AveragePooling2D(pool_size=(2, 2), strides=(2, 2))(x)\n    \n    # 64\n    x = Conv2D(32, kernel_size=(3, 3), strides=(1, 1), padding='same')(x)\n    for _ in range(3):\n        x = res_block(x, 32)\n    x = AveragePooling2D(pool_size=(2, 2), strides=(2, 2))(x)\n    \n    # 32\n    x = Conv2D(64, kernel_size=(3, 3), strides=(1, 1), padding='same')(x)\n    for _ in range(3):\n        x = res_block(x, 64)\n    x = AveragePooling2D(pool_size=(2, 2), strides=(2, 2))(x)    \n    \n    # 16\n    x = Conv2D(128, kernel_size=(3, 3), strides=(1, 1), padding='same')(x)\n    for _ in range(3):\n        x = res_block(x, 128)\n        \n    # 16\n    x = GlobalAveragePooling2D()(x)\n    \n    inp2 = Input(shape=(4,))\n    x2 = tf.keras.layers.GaussianNoise(0.2)(inp2)\n    x = Concatenate()([x, x2]) \n    x = Dropout(0.6)(x) \n    x = Dense(1)(x)\n    #x2 = Dense(1)(x)\n    return Model([inp, inp2] , x)","metadata":{"execution":{"iopub.status.busy":"2021-11-27T06:27:47.592831Z","iopub.execute_input":"2021-11-27T06:27:47.593501Z","iopub.status.idle":"2021-11-27T06:27:47.724647Z","shell.execute_reply.started":"2021-11-27T06:27:47.59346Z","shell.execute_reply":"2021-11-27T06:27:47.723723Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"model = get_model() \nmodel.summary() ","metadata":{"execution":{"iopub.status.busy":"2021-11-27T06:27:47.726082Z","iopub.execute_input":"2021-11-27T06:27:47.726486Z","iopub.status.idle":"2021-11-27T06:27:50.90164Z","shell.execute_reply.started":"2021-11-27T06:27:47.726415Z","shell.execute_reply":"2021-11-27T06:27:50.900905Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from tensorflow_addons.optimizers import RectifiedAdam\n\nmodel.compile(optimizer=tf.keras.optimizers.Adam(learning_rate=0.001), loss='mae') ","metadata":{"execution":{"iopub.status.busy":"2021-11-27T06:27:50.905315Z","iopub.execute_input":"2021-11-27T06:27:50.905582Z","iopub.status.idle":"2021-11-27T06:27:51.031404Z","shell.execute_reply.started":"2021-11-27T06:27:50.905556Z","shell.execute_reply":"2021-11-27T06:27:51.030696Z"},"trusted":true},"execution_count":null,"outputs":[]},{"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) ","metadata":{"execution":{"iopub.status.busy":"2021-11-27T06:27:51.032615Z","iopub.execute_input":"2021-11-27T06:27:51.032938Z","iopub.status.idle":"2021-11-27T06:27:51.586436Z","shell.execute_reply.started":"2021-11-27T06:27:51.032901Z","shell.execute_reply":"2021-11-27T06:27:51.585558Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import seaborn as sns\n\nsns.distplot(list(A.values()));","metadata":{"execution":{"iopub.status.busy":"2021-11-27T06:27:51.588691Z","iopub.execute_input":"2021-11-27T06:27:51.588962Z","iopub.status.idle":"2021-11-27T06:27:51.84361Z","shell.execute_reply.started":"2021-11-27T06:27:51.588937Z","shell.execute_reply":"2021-11-27T06:27:51.842664Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"er = tf.keras.callbacks.EarlyStopping(\n    monitor=\"val_loss\",\n    min_delta=1e-3,\n    patience=5,\n    verbose=0,\n    mode=\"auto\",\n    baseline=None,\n    restore_best_weights=True,\n)","metadata":{"execution":{"iopub.status.busy":"2021-11-27T06:27:51.845215Z","iopub.execute_input":"2021-11-27T06:27:51.845814Z","iopub.status.idle":"2021-11-27T06:27:51.851707Z","shell.execute_reply.started":"2021-11-27T06:27:51.845774Z","shell.execute_reply":"2021-11-27T06:27:51.850842Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"model.fit_generator(IGenerator(keys=tr_p, \n                               a = A, \n                               tab = TAB), \n                    steps_per_epoch = 200,\n                    validation_data=IGenerator(keys=vl_p, \n                               a = A, \n                               tab = TAB),\n                    validation_steps = 20, \n                    callbacks = [er], \n                    epochs=15)","metadata":{"execution":{"iopub.status.busy":"2021-11-27T06:27:51.853472Z","iopub.execute_input":"2021-11-27T06:27:51.853917Z","iopub.status.idle":"2021-11-27T06:38:37.605462Z","shell.execute_reply.started":"2021-11-27T06:27:51.853877Z","shell.execute_reply":"2021-11-27T06:38:37.6047Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def score(fvc_true, fvc_pred, sigma):\n    sigma_clip = np.maximum(sigma, 70)\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)","metadata":{"execution":{"iopub.status.busy":"2021-11-27T06:38:37.606987Z","iopub.execute_input":"2021-11-27T06:38:37.60832Z","iopub.status.idle":"2021-11-27T06:38:37.615278Z","shell.execute_reply.started":"2021-11-27T06:38:37.608283Z","shell.execute_reply":"2021-11-27T06:38:37.614462Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from tqdm.notebook import tqdm\n\nmetric = []\nfor 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        for i in os.listdir(f'../input/osic-pulmonary-fibrosis-progression/train/{p}/'):\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        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))","metadata":{"execution":{"iopub.status.busy":"2021-11-27T06:38:37.617007Z","iopub.execute_input":"2021-11-27T06:38:37.617503Z","iopub.status.idle":"2021-11-27T06:47:15.757134Z","shell.execute_reply.started":"2021-11-27T06:38:37.617465Z","shell.execute_reply":"2021-11-27T06:47:15.756258Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Predict","metadata":{}},{"cell_type":"code","source":"q = (np.argmin(metric) + 1)/ 10\nq","metadata":{"execution":{"iopub.status.busy":"2021-11-27T06:47:15.75882Z","iopub.execute_input":"2021-11-27T06:47:15.759423Z","iopub.status.idle":"2021-11-27T06:47:15.767598Z","shell.execute_reply.started":"2021-11-27T06:47:15.759383Z","shell.execute_reply":"2021-11-27T06:47:15.766683Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sub = pd.read_csv('../input/osic-pulmonary-fibrosis-progression/sample_submission.csv') \nsub.head() ","metadata":{"execution":{"iopub.status.busy":"2021-11-27T06:47:15.769136Z","iopub.execute_input":"2021-11-27T06:47:15.7698Z","iopub.status.idle":"2021-11-27T06:47:15.817914Z","shell.execute_reply.started":"2021-11-27T06:47:15.76976Z","shell.execute_reply":"2021-11-27T06:47:15.817279Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"test = pd.read_csv('../input/osic-pulmonary-fibrosis-progression/test.csv') \ntest.head()\n#len(test)","metadata":{"execution":{"iopub.status.busy":"2021-11-27T06:59:15.782119Z","iopub.execute_input":"2021-11-27T06:59:15.782463Z","iopub.status.idle":"2021-11-27T06:59:15.800092Z","shell.execute_reply.started":"2021-11-27T06:59:15.782427Z","shell.execute_reply":"2021-11-27T06:59:15.79922Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"A_test, B_test, P_test,W, FVC= {}, {}, {},{},{} \nSTD, WEEK = {}, {} \nfor p in test.Patient.unique():\n    x = [] \n    tab = [] \n    for i in os.listdir(f'../input/osic-pulmonary-fibrosis-progression/test/{p}/'):\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    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]","metadata":{"execution":{"iopub.status.busy":"2021-11-27T06:47:15.845116Z","iopub.execute_input":"2021-11-27T06:47:15.845678Z","iopub.status.idle":"2021-11-27T06:47:41.432288Z","shell.execute_reply.started":"2021-11-27T06:47:15.845639Z","shell.execute_reply":"2021-11-27T06:47:41.431297Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"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    ","metadata":{"execution":{"iopub.status.busy":"2021-11-27T06:47:41.433657Z","iopub.execute_input":"2021-11-27T06:47:41.433993Z","iopub.status.idle":"2021-11-27T06:47:44.224635Z","shell.execute_reply.started":"2021-11-27T06:47:41.433958Z","shell.execute_reply":"2021-11-27T06:47:44.223866Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sub.head()","metadata":{"execution":{"iopub.status.busy":"2021-11-27T06:47:44.225867Z","iopub.execute_input":"2021-11-27T06:47:44.226208Z","iopub.status.idle":"2021-11-27T06:47:44.238789Z","shell.execute_reply.started":"2021-11-27T06:47:44.226163Z","shell.execute_reply":"2021-11-27T06:47:44.237717Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sub[[\"Patient_Week\",\"FVC\",\"Confidence\"]].to_csv(\"submission.csv\", index=False)","metadata":{"execution":{"iopub.status.busy":"2021-11-27T06:47:44.240792Z","iopub.execute_input":"2021-11-27T06:47:44.241591Z","iopub.status.idle":"2021-11-27T06:47:44.642198Z","shell.execute_reply.started":"2021-11-27T06:47:44.241552Z","shell.execute_reply":"2021-11-27T06:47:44.641461Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"submi = pd.read_csv('submission.csv') \nsubmi.head()","metadata":{"execution":{"iopub.status.busy":"2021-11-27T06:47:44.644253Z","iopub.execute_input":"2021-11-27T06:47:44.644865Z","iopub.status.idle":"2021-11-27T06:47:44.661987Z","shell.execute_reply.started":"2021-11-27T06:47:44.644825Z","shell.execute_reply":"2021-11-27T06:47:44.661298Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{"trusted":true},"execution_count":null,"outputs":[]}]}