{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.10.14","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"gpu","dataSources":[{"sourceId":71549,"databundleVersionId":8561470,"sourceType":"competition"},{"sourceId":9459407,"sourceType":"datasetVersion","datasetId":5750674}],"dockerImageVersionId":30761,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":true}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"import pydicom as dicom\nimport matplotlib.pylab as plt\nimport pandas as pd\nimport numpy as np\nimport os\nimport torch\nimport torch.nn as nn\nimport torch.optim as optim\nimport torch.nn.functional as F\nfrom torch.utils.data import DataLoader, Dataset, WeightedRandomSampler\nimport random\nfrom sklearn.model_selection import StratifiedShuffleSplit\nfrom concurrent.futures import ThreadPoolExecutor\nfrom concurrent.futures import as_completed\nimport gzip\nimport re\nimport torchvision.transforms.functional as FF","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:25:18.852086Z","iopub.execute_input":"2024-10-06T12:25:18.852782Z","iopub.status.idle":"2024-10-06T12:25:24.507366Z","shell.execute_reply.started":"2024-10-06T12:25:18.852742Z","shell.execute_reply":"2024-10-06T12:25:24.506463Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\"\"\" \nDefining some paths to make it easier to change the file structure.\n\"\"\"\n\ndata_folder = \"../input/rsna-serial/rsna-2024/\"\ninput_folder = \"../input/rsna-2024-lumbar-spine-degenerative-classification/\"\nimage_folder = \"train_images\"\ntensors_folder = \"serial_train\"\n\ndf_s_st1_path = os.path.join(data_folder,\"df_s_st1.csv\")\ndf_s_st2_path = os.path.join(data_folder,\"df_s_st2.csv\")\ndf_s_at2_path = os.path.join(data_folder,\"df_s_at2.csv\")\ndf_desc = pd.read_csv(os.path.join(data_folder,\"train_series_descriptions.csv\"))\n\n# Needed only to serialize\n# This was already done into the dataset \"rsna-serial\"\ndf_st1_path = os.path.join(\"df_st1.csv\")\ndf_st2_path = os.path.join(\"df_st2.csv\")\ndf_at2_path = os.path.join(\"df_at2.csv\")","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:25:24.509069Z","iopub.execute_input":"2024-10-06T12:25:24.509576Z","iopub.status.idle":"2024-10-06T12:25:24.541378Z","shell.execute_reply.started":"2024-10-06T12:25:24.509542Z","shell.execute_reply":"2024-10-06T12:25:24.540306Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\"\"\"\nHere I am separating the descriptions into each imagery type\n\"\"\"\ndf_st1 = df_desc[df_desc[\"series_description\"] == \"Sagittal T1\"]\ndf_st2 = df_desc[df_desc[\"series_description\"] == \"Sagittal T2/STIR\"]\ndf_at2 = df_desc[df_desc[\"series_description\"] == \"Axial T2\"]\n\nprint(\"S_T1 series number: \",len(df_at2))\nprint(\"S_T2 series number: \",len(df_st1))\nprint(\"A_T1 series number: \",len(df_st2))\n","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:25:24.542737Z","iopub.execute_input":"2024-10-06T12:25:24.543177Z","iopub.status.idle":"2024-10-06T12:25:24.563703Z","shell.execute_reply.started":"2024-10-06T12:25:24.543133Z","shell.execute_reply":"2024-10-06T12:25:24.562693Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\"\"\"\nSome functions\n\"\"\"\ndef th_get_size(image_path,study_id,series_id,img_name):\n    img_data = {}\n    ds = dicom.dcmread(image_path, stop_before_pixels=True)\n    img_data[\"study_id\"] = study_id\n    img_data[\"series_id\"] = series_id\n    img_data[\"img_name\"] = img_name\n    img_data[\"cols\"] = ds.Rows\n    img_data[\"rows\"] = ds.Columns\n    img_data[\"max\"] = np.max([ds.Rows,ds.Columns])\n    return img_data\n\ndef get_images_number(df):\n    img_number = []\n    for x in df.iterrows():\n        serie_path = os.path.join(data_folder,image_folder,str(x[1].iloc[0]),str(x[1].iloc[1]))\n        images = os.listdir(serie_path)\n        img_number.append(len(images))\n    df[\"image_number\"] = img_number\n\ndef get_sizes_df(df,num_threads=4):\n    data = {\"study_id\":[],\"series_id\":[],\"img_name\":[],\"cols\":[],\"rows\":[],\"max\":[]}\n    futures = []\n    img_number = []\n    with ThreadPoolExecutor(num_threads) as executor:\n        for x in df.iterrows():\n            serie_path = os.path.join(data_folder,image_folder,str(x[1].iloc[0]),str(x[1].iloc[1]))\n            images = os.listdir(serie_path)\n            img_number.append(len(images))\n            for img_name in images:\n                image_path = os.path.join(serie_path,img_name)\n                futures.append(executor.submit(th_get_size,image_path,x[1].iloc[0],x[1].iloc[1],img_name))\n    for future in as_completed(futures):\n        dick = future.result()\n        data[\"study_id\"].append(dick[\"study_id\"])\n        data[\"series_id\"].append(dick[\"series_id\"]) \n        data[\"img_name\"].append(dick[\"img_name\"])\n        data[\"cols\"].append(dick[\"cols\"])\n        data[\"rows\"].append(dick[\"rows\"])\n        data[\"max\"].append(dick[\"max\"])\n\n    return pd.DataFrame.from_dict(data,orient=\"columns\")\n\ndef plot_image(df,img_id):\n    study_id = df.loc[img_id,\"study_id\"]\n    series_id = df.loc[img_id,\"series_id\"]\n    img_name = df.loc[img_id,\"img_name\"]\n    image_path = os.path.join(input_folder,image_folder,str(study_id),str(series_id),str(img_name))\n    ds = dicom.dcmread(image_path)\n    plt.imshow(ds.pixel_array)\n    plt.title(f\"{ds.pixel_array.shape}\")\n    df = pd.read_csv(os.path.join(data_folder,\"train_label_coordinates.csv\"))\n    df = df[(df[\"study_id\"] == study_id) & (df[\"series_id\"] == series_id)]\n    df_coords = df[df[\"instance_number\"] == int(img_name.split(\".\")[0])].reset_index(drop=True)\n    if len(df_coords) != 0:\n        plt.scatter(df_coords[\"x\"],df_coords[\"y\"],s=50,facecolors='none', edgecolors='r',marker=\"o\")\n\ndef plot_sizes_distributions(df_arr):\n    fig,ax = plt.subplots(1,3,figsize=(10,3))\n    fig.suptitle(\"Sizes distribution for each imagery type\")\n    for i,df_obj in enumerate(df_arr):\n        df_dist = df_obj[\"df\"][\"max\"].value_counts().reset_index().sort_values(by=\"max\")\n        ax[i].bar(df_dist[\"max\"],df_dist[\"count\"],width=20)\n        ax[i].set_xlabel(\"Max edge size\")\n        ax[i].set_ylabel(\"Counts\")\n        ax[i].set_title(df_obj[\"name\"]+\"\\n\"+f\"Min:{df_obj['df']['max'].min()}\"+\"\\n\"+f\"Max:{df_obj['df']['max'].max()}\")\n    plt.tight_layout()","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:25:24.566737Z","iopub.execute_input":"2024-10-06T12:25:24.567301Z","iopub.status.idle":"2024-10-06T12:25:24.588950Z","shell.execute_reply.started":"2024-10-06T12:25:24.567262Z","shell.execute_reply":"2024-10-06T12:25:24.587931Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\"\"\"\nHere I am getting each images sizes and number of images for each imagery type\n\"\"\"\nnum_threads = 4\n\n\nif os.path.exists(df_s_st1_path):\n    df_s_st1 = pd.read_csv(df_s_st1_path)\nelse:\n    df_s_st1 = get_sizes_df(df_st1,num_threads)\n    df_s_st1.to_csv(df_s_st1_path)\n\nif os.path.exists(df_s_st2_path):\n    df_s_st2 = pd.read_csv(df_s_st2_path)\nelse:\n    df_s_st2 = get_sizes_df(df_st2,num_threads)\n    df_s_st2.to_csv(df_s_st2_path)\n\nif os.path.exists(df_s_at2_path):\n    df_s_at2 = pd.read_csv(df_s_at2_path)\nelse:\n    df_s_at2 = get_sizes_df(df_at2,num_threads)\n    df_s_at2.to_csv(df_s_at2_path)\n    \n\n","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:25:24.590138Z","iopub.execute_input":"2024-10-06T12:25:24.590480Z","iopub.status.idle":"2024-10-06T12:25:24.825063Z","shell.execute_reply.started":"2024-10-06T12:25:24.590447Z","shell.execute_reply":"2024-10-06T12:25:24.824066Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\"\"\"\nHere we can see the percentage of images that aren't perfect squares \n\"\"\"\ndf_arr = [\n    {\"name\":\"Saggital T1\",\"df\":df_s_st1},\n    {\"name\":\"Saggital T2/STIR\",\"df\":df_s_st2},\n    {\"name\":\"Axial T2\",\"df\":df_s_at2}\n]\nfor df in df_arr:\n    not_square = df[\"df\"][df[\"df\"][\"cols\"] != df[\"df\"][\"rows\"]]\n    ratio = len(not_square)/len(df_s_st1)\n    print(df[\"name\"], f\"({len(df['df'])} images)\")\n    print(\"Percentage of non-square:\",ratio,\"\\n\")\n","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:25:24.826574Z","iopub.execute_input":"2024-10-06T12:25:24.827001Z","iopub.status.idle":"2024-10-06T12:25:24.838280Z","shell.execute_reply.started":"2024-10-06T12:25:24.826961Z","shell.execute_reply":"2024-10-06T12:25:24.837193Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_plot = df_s_st1[(df_s_st1[\"cols\"]!=df_s_st1[\"rows\"])&\n         (np.abs(df_s_st1[\"cols\"]-df_s_st1[\"rows\"])>100)].sort_values(by=[\"series_id\",\"img_name\"]).reset_index(drop=True)\n\nplot_image(df_plot,2)","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:25:24.839584Z","iopub.execute_input":"2024-10-06T12:25:24.839936Z","iopub.status.idle":"2024-10-06T12:25:27.757564Z","shell.execute_reply.started":"2024-10-06T12:25:24.839874Z","shell.execute_reply":"2024-10-06T12:25:27.756554Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plot_image(df_plot,68)","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:25:27.758688Z","iopub.execute_input":"2024-10-06T12:25:27.758977Z","iopub.status.idle":"2024-10-06T12:25:28.297811Z","shell.execute_reply.started":"2024-10-06T12:25:27.758947Z","shell.execute_reply":"2024-10-06T12:25:28.296851Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plot_image(df_plot,171)","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:25:28.298967Z","iopub.execute_input":"2024-10-06T12:25:28.299248Z","iopub.status.idle":"2024-10-06T12:25:28.668159Z","shell.execute_reply.started":"2024-10-06T12:25:28.299217Z","shell.execute_reply":"2024-10-06T12:25:28.667224Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_plot = df_s_at2[(df_s_at2[\"cols\"]!=df_s_at2[\"rows\"])&\n         (np.abs(df_s_at2[\"cols\"]-df_s_at2[\"rows\"])>50)].sort_values(by=[\"series_id\",\"img_name\"]).reset_index(drop=True)\nplot_image(df_plot,8)","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:25:28.672798Z","iopub.execute_input":"2024-10-06T12:25:28.673298Z","iopub.status.idle":"2024-10-06T12:25:29.182857Z","shell.execute_reply.started":"2024-10-06T12:25:28.673253Z","shell.execute_reply":"2024-10-06T12:25:29.181936Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plot_image(df_plot,33)","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:25:29.184222Z","iopub.execute_input":"2024-10-06T12:25:29.184866Z","iopub.status.idle":"2024-10-06T12:25:29.597758Z","shell.execute_reply.started":"2024-10-06T12:25:29.184821Z","shell.execute_reply":"2024-10-06T12:25:29.596826Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plot_image(df_plot,77)","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:25:29.599103Z","iopub.execute_input":"2024-10-06T12:25:29.599519Z","iopub.status.idle":"2024-10-06T12:25:30.046123Z","shell.execute_reply.started":"2024-10-06T12:25:29.599464Z","shell.execute_reply":"2024-10-06T12:25:30.045095Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plot_image(df_plot,103)","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:25:30.047365Z","iopub.execute_input":"2024-10-06T12:25:30.047670Z","iopub.status.idle":"2024-10-06T12:25:30.554430Z","shell.execute_reply.started":"2024-10-06T12:25:30.047637Z","shell.execute_reply":"2024-10-06T12:25:30.553464Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plot_image(df_plot,339)","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:25:30.555765Z","iopub.execute_input":"2024-10-06T12:25:30.556204Z","iopub.status.idle":"2024-10-06T12:25:31.128695Z","shell.execute_reply.started":"2024-10-06T12:25:30.556130Z","shell.execute_reply":"2024-10-06T12:25:31.127684Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plot_sizes_distributions(df_arr)","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:25:31.130045Z","iopub.execute_input":"2024-10-06T12:25:31.130369Z","iopub.status.idle":"2024-10-06T12:25:32.282760Z","shell.execute_reply.started":"2024-10-06T12:25:31.130333Z","shell.execute_reply":"2024-10-06T12:25:32.281772Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\"\"\"\nSerialize data into tensors\n\"\"\"\ndef extract_number(filename):\n    match = re.search(r'(\\d+)', filename)\n    return int(match.group(0)) if match else float('inf')\n\ndef serialize(df_d,df_s,imagery_type):\n    resize_size = round(df_s[\"max\"].mean())\n    resize_depth = round(df_d[\"image_number\"].mean())\n    for i in range(0,len(df_d)):\n        df = df_d.iloc[i,:]\n        serie_path = os.path.join(data_folder,image_folder,str(df[\"study_id\"]),str(df[\"series_id\"]))\n        max_rows = df_s[df_s[\"series_id\"] == df[\"series_id\"]][\"rows\"].max()\n        max_cols = df_s[df_s[\"series_id\"] == df[\"series_id\"]][\"cols\"].max()\n        serie_tensor = []\n        images_sorted = sorted(os.listdir(serie_path), key=extract_number)\n        for image in images_sorted:\n            image_path = os.path.join(serie_path,image)\n            dicom_image = dicom.dcmread(image_path).pixel_array\n            dicom_image = torch.tensor(dicom_image,dtype=torch.int16)\n            if dicom_image.shape[0] < max_rows:\n                padding = (max_rows - dicom_image.shape[0])\n                padding_left = padding // 2\n                padding_right = padding - padding_left\n                dicom_image = F.pad(dicom_image,(0,0,padding_left,padding_right),value=0)\n            if dicom_image.shape[1] < max_cols:\n                padding = (max_cols - dicom_image.shape[1])\n                padding_top = padding // 2\n                padding_bottom = padding - padding_top\n                dicom_image = F.pad(dicom_image,(padding_top,padding_bottom,0,0),value=0)\n            dicom_image = dicom_image.unsqueeze(0)\n            dicom_image = FF.resize(dicom_image,resize_size)\n            serie_tensor.append(dicom_image)\n\n        stack = torch.stack(serie_tensor,dim=1)\n        stack_resized = []\n        for i in range(stack.shape[2]):\n            stack_resized.append(FF.resize(stack[0,:,i,:].unsqueeze(0),size=(resize_depth,resize_size)))\n        stack_resized = torch.stack(stack_resized,dim=2)    \n        study_dir = os.path.join(tensors_folder,imagery_type,str(df[\"study_id\"]))\n        if not os.path.exists(study_dir):\n            os.makedirs(study_dir)\n        path = os.path.join(study_dir,str(df[\"series_id\"])+\".pth.gz\")\n        with gzip.open(path, 'wb') as f:\n            torch.save(stack_resized,f)\n\n# To launch in order to serialize\n\"\"\"\nif os.path.exists(df_st1_path):\n    df_st1 = pd.read_csv(df_st1_path)\nelse:\n    get_images_number(df_st1)\n    df_st1.to_csv(df_st1_path)\n\nif os.path.exists(df_st2_path):\n    df_st2 = pd.read_csv(df_st2_path)\nelse:\n    get_images_number(df_st2)\n    df_st2.to_csv(df_st2_path)\n\nif os.path.exists(df_at2_path):\n    df_at2 = pd.read_csv(df_at2_path)\nelse:\n    get_images_number(df_at2)\n    df_at2.to_csv(df_at2_path)\n\nserialize(df_st1,df_s_st1,\"st1\")\nserialize(df_st2,df_s_st2,\"st2\")\nserialize(df_at2,df_s_at2,\"at2\")\n\"\"\"","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:25:32.284649Z","iopub.execute_input":"2024-10-06T12:25:32.285275Z","iopub.status.idle":"2024-10-06T12:25:32.306043Z","shell.execute_reply.started":"2024-10-06T12:25:32.285225Z","shell.execute_reply":"2024-10-06T12:25:32.305134Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\"\"\"We can see some vertebraes, hopefully the 3D CNN will get some extra informations along this dimension\"\"\"\npath =os.path.join(data_folder,tensors_folder,\"st1/4003253/1054713880.pth\")\ntensor = torch.load(path,weights_only=False)\n\ntranche = tensor[0,:,:,200]\nprint(tranche.shape)\nplt.imshow(tranche.T.numpy())","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:25:32.307367Z","iopub.execute_input":"2024-10-06T12:25:32.307674Z","iopub.status.idle":"2024-10-06T12:25:32.646460Z","shell.execute_reply.started":"2024-10-06T12:25:32.307643Z","shell.execute_reply":"2024-10-06T12:25:32.645571Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\"\"\"\npreparing to group using df_labels\n\"\"\"\ndf_labels = pd.read_csv(os.path.join(data_folder,\"train.csv\"))\nonehot_map = {\n    \"Normal/Mild\":0,\n    \"Moderate\":1,\n    \"Severe\":2,\n    \"NaN\":3\n}\n\ncolabel = df_labels.columns.difference([\"study_id\"])\ndf_labels[\"colabel\"] = [0.0 for x in df_labels.iterrows()]\nfor col in colabel:\n    df_labels[col] = df_labels[col].fillna(\"NaN\").map(onehot_map)\n    df_labels[\"colabel\"] += df_labels[col]\n\n# onehotting\ncolabel = df_labels[\"colabel\"]\nstudy_id = df_labels[\"study_id\"]\ndf_labels = df_labels.drop([\"colabel\",\"study_id\"],axis=1)\nfor col in df_labels.columns:\n    cols = {}\n    for i in range(0,4):\n        cols[col+\"_\"+str(i)] = df_labels[col].map(lambda x: 1 if x == i else 0)\n    new_cols = pd.DataFrame(cols)\n    df_labels = pd.concat([df_labels,new_cols],axis=1) \n    df_labels = df_labels.drop(col,axis=1)\ndf_labels = df_labels.astype(float)\ndf_labels[\"study_id\"] = study_id\ndf_labels[\"colabel\"] = colabel\ndf_counts = df_labels[\"colabel\"].value_counts().reset_index()\n\n# grouping colabels with only one element (usefull to use a groupsplit strategy per study)\ndf_mins = df_counts[df_counts[\"count\"]==1]\ndf_counts[\"g_colabel\"] = df_counts[\"colabel\"]\ndf_counts[\"g_colabel\"] = [df_mins.loc[df_mins[\"colabel\"] == df_mins[\"colabel\"].min(),\"colabel\"].item() if x[1].iloc[1] == 1 else x[1].iloc[0] for x in df_counts.iterrows() ]\nmap_colabel = df_counts.drop(\"count\",axis=1)\nmap_colabel = {x[1].iloc[0]:x[1].iloc[1] for x in map_colabel.iterrows()}\ndf_labels[\"colabel\"] = df_labels[\"colabel\"].map(map_colabel)\n\n# creating a list of classes distribution\ncriterion_weights = []\ncols = df_labels.columns.difference([\"study_id\",\"colabel\"],sort=False)\nfor i in range(0,25):\n    target_cols = cols[i*4:i*4+4]\n    df_distrib = df_labels[target_cols].value_counts().reset_index()\n    summ = df_distrib[\"count\"].sum()\n    for i,col in enumerate(target_cols):\n        ratio = df_distrib[df_distrib[col] == 1.0][\"count\"].item()/summ\n        if i == 3:\n            final_ratio = 0.99\n        else:\n            final_ratio = ratio \n        criterion_weights.append(final_ratio)\n        \ncriterion_weights = 1/torch.tensor(criterion_weights)","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:25:32.647906Z","iopub.execute_input":"2024-10-06T12:25:32.648575Z","iopub.status.idle":"2024-10-06T12:25:33.158732Z","shell.execute_reply.started":"2024-10-06T12:25:32.648528Z","shell.execute_reply":"2024-10-06T12:25:33.157932Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\"\"\"\nFunction to distribute targets equally across splits\n\"\"\"\ndef get_splits(df_labels,df_d,n_splits=5,random_state=42):\n    sss = StratifiedShuffleSplit(n_splits=n_splits,test_size=0.2,random_state=random_state)\n    splits = []\n    present_studies = list(pd.unique(df_d[\"study_id\"]))\n    df_l = df_labels[df_labels[\"study_id\"].isin(present_studies)].reset_index(drop=True)\n    for train,val in sss.split(df_l,df_l[\"colabel\"]):\n        splits.append({\"train\":df_l.loc[train,\"study_id\"].to_list(),\"val\":df_l.loc[val,\"study_id\"].to_list()})\n    return splits","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:25:33.159830Z","iopub.execute_input":"2024-10-06T12:25:33.160144Z","iopub.status.idle":"2024-10-06T12:25:33.167240Z","shell.execute_reply.started":"2024-10-06T12:25:33.160111Z","shell.execute_reply":"2024-10-06T12:25:33.166045Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\"\"\"\nclass to build the train and validation datasets\n\"\"\"\nclass DataGetter(Dataset):\n    def __init__(self, df_labels, df_transforms = None, imagery=\"st1\"):\n        self.imagery = imagery\n        self.df_labels = df_labels\n        self.df_transforms = df_transforms\n        self.series_paths = []\n        if type(self.df_transforms) == type(pd.DataFrame()):\n            self.training = True\n        else:\n            self.training = False\n\n    def __len__(self):\n        return len(self.series_paths)\n    \n    def __getitem__(self, idx):\n        serie_path = self.series_paths[idx]\n        tensor = torch.load(serie_path,weights_only=False)\n        tensor = tensor.to(torch.float)\n        mint = torch.min(tensor)\n        tensor = (tensor-mint)/(torch.max(tensor)-mint)\n        if self.training:\n            mean = self.df_transforms[self.df_transforms[\"imagery\"] == self.imagery][\"mean\"]\n            std = self.df_transforms[self.df_transforms[\"imagery\"] == self.imagery][\"std\"]\n            tensor = ((tensor-torch.tensor(mean))/torch.tensor(std)).to(torch.float32)\n        path_split = os.path.normpath(serie_path).split(os.path.sep)\n        labels = self.df_labels[self.df_labels[\"study_id\"] == int(path_split[-2])].drop([\"study_id\",\"colabel\"],axis=1)\n        labels = torch.tensor(labels.reset_index(drop=True).to_numpy())\n        return tensor,labels\n\n    def load_series_paths(self,study_ids):\n        self.series_paths = []\n        base_path = os.path.join(data_folder,tensors_folder,self.imagery)\n        for study in study_ids:\n            study = str(study)\n            study_path = os.path.join(base_path, study)\n            for serie in os.listdir(study_path):\n                self.series_paths.append(os.path.join(study_path, serie))","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:25:33.168747Z","iopub.execute_input":"2024-10-06T12:25:33.169222Z","iopub.status.idle":"2024-10-06T12:25:33.181837Z","shell.execute_reply.started":"2024-10-06T12:25:33.169179Z","shell.execute_reply":"2024-10-06T12:25:33.180912Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\"\"\"\nCreating splits and\nCalculating mean and stddev for training data\n\"\"\"\nif not os.path.exists(os.path.join(data_folder,\"df_transforms.csv\")):\n    splits_st1 = get_splits(df_labels,df_st1)\n    data_st1 = DataGetter(df_labels,\"st1\")\n    data_st1.load_series_paths(splits_st1[0][\"train\"])\n    loader_st1 = torch.utils.data.DataLoader(data_st1, batch_size=1, shuffle=False, num_workers=0)\n\n    splits_st2 = get_splits(df_labels,df_st2)\n    data_st2 = DataGetter(df_labels,\"st2\")\n    data_st2.load_series_paths(splits_st2[0][\"train\"])\n    loader_st2 = torch.utils.data.DataLoader(data_st2, batch_size=1, shuffle=False, num_workers=0)\n\n    splits_at2 = get_splits(df_labels,df_at2)\n    data_at2 = DataGetter(df_labels,\"at2\")\n    data_at2.load_series_paths(splits_at2[0][\"train\"])\n    loader_at2 = torch.utils.data.DataLoader(data_at2, batch_size=1, shuffle=False, num_workers=0)\n\n    data = {\"imagery\":[],\"mean\":[],\"std\":[]}\n    mean = 0\n    std = 0\n    for img,_ in loader_st1:\n        mean += img.mean()\n        std += img.std()\n        print(img.shape,mean,std)\n    mean /= len(data_st1)\n    std /= len(data_st1)\n    data[\"imagery\"].append(\"st1\")\n    data[\"mean\"].append(mean.item())\n    data[\"std\"].append(std.item())\n\n    mean = 0\n    std = 0\n    index = 0\n    for img,_ in loader_st2:\n        mean += img.mean()\n        std += img.std()\n    mean /= len(data_st2)\n    std /= len(data_st2)\n    data[\"imagery\"].append(\"st2\")\n    data[\"mean\"].append(mean.item())\n    data[\"std\"].append(std.item())\n\n    mean = 0\n    std = 0\n    index = 0\n    for img,_ in loader_at2:\n        mean += img.mean()\n        std += img.std()\n    mean /= len(data_at2)\n    std /= len(data_at2)\n    data[\"imagery\"].append(\"at2\")\n    data[\"mean\"].append(mean.item())\n    data[\"std\"].append(std.item())\n    \n    d_split_st1 = {\"study_id\":[],\"split\":[]}\n    for study in splits_st1[0][\"train\"]:\n        d_split_st1[\"study_id\"].append(study)\n        d_split_st1[\"split\"].append(\"train\")\n    for study in splits_st1[0][\"val\"]:\n        d_split_st1[\"study_id\"].append(study)\n        d_split_st1[\"split\"].append(\"val\")\n\n    d_split_st2 = {\"study_id\":[],\"split\":[]}\n    for study in splits_st2[0][\"train\"]:\n        d_split_st2[\"study_id\"].append(study)\n        d_split_st2[\"split\"].append(\"train\")\n    for study in splits_st2[0][\"val\"]:\n        d_split_st2[\"study_id\"].append(study)\n        d_split_st2[\"split\"].append(\"val\")\n\n    d_split_at2 = {\"study_id\":[],\"split\":[]}\n    for study in splits_at2[0][\"train\"]:\n        d_split_at2[\"study_id\"].append(study)\n        d_split_at2[\"split\"].append(\"train\")\n    for study in splits_at2[0][\"val\"]:\n        d_split_at2[\"study_id\"].append(study)\n        d_split_at2[\"split\"].append(\"val\")\n\n    df_split_st1 = pd.DataFrame(data=d_split_st1)\n    df_split_st1.to_csv(os.path.join(data_folder,\"df_split_st1.csv\"),index=False)\n    df_split_st2 = pd.DataFrame(data=d_split_st2)\n    df_split_st2.to_csv(os.path.join(data_folder,\"df_split_st2.csv\"),index=False)\n    df_split_at2 = pd.DataFrame(data=d_split_at2)\n    df_split_at2.to_csv(os.path.join(data_folder,\"df_split_at2.csv\"),index=False)\n    df_transforms = pd.DataFrame(data=data)\n    df_transforms.to_csv(os.path.join(data_folder,\"df_transforms.csv\"),index=False)\nelse:\n    df_transforms = pd.read_csv(os.path.join(data_folder,\"df_transforms.csv\"))\n    df_split_st1 = pd.read_csv(os.path.join(data_folder,\"df_split_st1.csv\"))\n    df_split_st2 = pd.read_csv(os.path.join(data_folder,\"df_split_st2.csv\"))\n    df_split_at2 = pd.read_csv(os.path.join(data_folder,\"df_split_at2.csv\"))","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:25:33.183269Z","iopub.execute_input":"2024-10-06T12:25:33.183576Z","iopub.status.idle":"2024-10-06T12:25:33.238254Z","shell.execute_reply.started":"2024-10-06T12:25:33.183543Z","shell.execute_reply":"2024-10-06T12:25:33.237450Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def calc(inpt,pad,kern,stride):\n    return (inpt+2*pad-(kern-1)-1)/stride +1\n    \ncalc(123,0,4,1)","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:25:33.239420Z","iopub.execute_input":"2024-10-06T12:25:33.239780Z","iopub.status.idle":"2024-10-06T12:25:33.247170Z","shell.execute_reply.started":"2024-10-06T12:25:33.239743Z","shell.execute_reply":"2024-10-06T12:25:33.246238Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"class ST1CNN(nn.Module):\n    def __init__(self):\n        super(ST1CNN, self).__init__()\n        self.actconv = nn.Tanh()\n        # 13 437 341\n        # Convolution layer 1\n        self.conv1 = nn.Conv3d(1, 2, kernel_size=(1,4,4), stride=(1,1,1), padding=(0,0,0))\n        # 13 434 338\n        self.pool1 = nn.AdaptiveAvgPool3d((13,320,250))\n        # 13 320 250\n        # actconv\n        # Convolution layer 2\n        self.conv2 = nn.Conv3d(2, 4, kernel_size=(1,4,4), stride=(1,1,1), padding=(0,0,0))\n        # 13 317 247\n        # actconv\n        self.pool2 = nn.AdaptiveAvgPool3d((13,234,182))\n        # 13 234 182\n        # Convolution layer 3\n        self.conv3 = nn.Conv3d(4, 8, kernel_size=(1,4,4), stride=(1,1,1), padding=(0,0,0))\n        # 13 231 179\n        # actconv\n        self.pool3 = nn.AdaptiveAvgPool3d((13,170,132))\n        # 13 170 132\n        # Convolution layer 4\n        self.conv4 = nn.Conv3d(8, 16, kernel_size=(1,4,4), stride=(1,1,1), padding=(0,0,0))\n        # 13 167 129\n        # actconv\n        self.pool4 = nn.AdaptiveAvgPool3d((13,123,95))\n        # 13 123 95\n        # Convolution layer 5\n        self.conv5 = nn.Conv3d(16, 32, kernel_size=(1,4,4), stride=(1,1,1), padding=(0,0,0))\n        # 13 120 92\n        # actconv\n        self.pool5 = nn.AdaptiveAvgPool3d((13,88,68))\n        # 13 88 68\n        # Convolution layer 6\n        self.conv6 = nn.Conv3d(32, 64, kernel_size=(1,3,3), stride=(1,1,1), padding=(0,0,0))\n        # 13 86 66\n        # actconv\n        self.pool6 = nn.AdaptiveAvgPool3d((13,64,48))\n        # 13 64 48\n        # Convolution layer 7\n        self.conv7 = nn.Conv3d(64, 64, kernel_size=(3,3,3), stride=(1,1,1), padding=(0,0,0))\n        # 11 62 46\n        # actconv\n        self.pool7 = nn.AdaptiveAvgPool3d((11,46,34))\n        # 11 46 34\n        # Convolution layer 8\n        self.conv8 = nn.Conv3d(64, 64, kernel_size=(3,3,3), stride=(1,1,1), padding=(0,0,0))\n        # 9 44 32\n        # actconv\n        self.pool8 = nn.AdaptiveAvgPool3d((9,32,24))\n        # 9 32 24\n        # Convolution layer 9\n        self.conv9 = nn.Conv3d(64, 64, kernel_size=(3,3,3), stride=(1,1,1), padding=(0,0,0))\n        # 7 30 22\n        # actconv\n        self.pool9 = nn.AdaptiveAvgPool3d((7,22,16))\n        # 7 22 16\n        # Convolution layer 10\n        self.conv10 = nn.Conv3d(64, 64, kernel_size=(3,3,3), stride=(1,1,1), padding=(0,0,0))\n        # 5 20 14\n        # actconv\n        self.pool10 = nn.AdaptiveAvgPool3d((5,14,10))\n        # 5 14 10\n        # Convolution layer 11\n        self.conv11 = nn.Conv3d(64, 64, kernel_size=(3,3,3), stride=(1,1,1), padding=(0,0,0))\n        # 3 12 8\n        # actconv\n        self.pool11 = nn.AdaptiveAvgPool3d((2,9,6))\n        # Fully connected layer\n        self.fc1 = nn.Linear(64*2*9*6, 100)\n        \n        #self.batchnorm = nn.BatchNorm3d(2, eps=1e-05, momentum=0.1, affine=True)\n        #self.maxpool = nn.MaxPool3d(kernel_size=(1,3,3),stride=(1,1,1),padding=(0,1,1))\n        #self.avgpool = nn.AvgPool3d(kernel_size=(1,3,3),stride=(1,1,1),padding=(0,1,1))\n        #self.attention_conv = nn.Conv3d(4,2,kernel_size=(1,3,3),stride=(1,1,1),padding=(0,1,1))\n\n    def forward(self, x):\n        # simulate negative padding\n        x = x[:,:,2:-2,16:-16,64:-64]\n        # Apply convolution 1\n        x = self.conv1(x)\n        x = self.pool1(x)\n        x = self.actconv(x)\n        # Apply convolution 2\n        x = self.conv2(x)\n        x = self.pool2(x)\n        x = self.actconv(x)\n        # Apply convolution 3\n        x = self.conv3(x)\n        x = self.pool3(x)\n        x = self.actconv(x)\n        # Apply convolution 4\n        x = self.conv4(x)\n        x = self.pool4(x)\n        x = self.actconv(x)\n        # Apply convolution 5\n        x = self.conv5(x)\n        x = self.pool5(x)\n        x = self.actconv(x)\n        # Apply convolution 6\n        x = self.conv6(x)\n        x = self.pool6(x)\n        x = self.actconv(x)\n        # Apply convolution 7\n        x = self.conv7(x)\n        x = self.pool7(x)\n        x = self.actconv(x)\n        # Apply convolution 8\n        x = self.conv8(x)\n        x = self.pool8(x)\n        x = self.actconv(x)\n        # Apply convolution 9\n        x = self.conv9(x)\n        x = self.pool9(x)\n        x = self.actconv(x)\n        # Apply convolution 10\n        x = self.conv10(x)\n        x = self.pool10(x)\n        x = self.actconv(x)\n        # Apply convolution 11\n        x = self.conv11(x)\n        x = self.pool11(x)\n        x = self.actconv(x)\n        x = torch.flatten(x,start_dim=1)\n        #x = self.batchnorm(x)\n        # Apply attention\n        #x1 = self.maxpool(x)\n        #x2 = self.avgpool(x)\n        #x = torch.cat([x1,x2],dim=1)\n        #x = self.attention_conv(x)\n        #x = self.actconv(x)\n        # Apply fully connected layers\n        return self.fc1(x)","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:32:25.484701Z","iopub.execute_input":"2024-10-06T12:32:25.485475Z","iopub.status.idle":"2024-10-06T12:32:25.513273Z","shell.execute_reply.started":"2024-10-06T12:32:25.485435Z","shell.execute_reply":"2024-10-06T12:32:25.512196Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\"\"\"\nFunction to ensure reproducibility\n\"\"\"\ndef set_random_seed(seed):\n    torch.manual_seed(seed)\n    np.random.seed(seed)\n    random.seed(seed)\n    \n    if torch.cuda.is_available():\n        torch.cuda.manual_seed(seed)\n        torch.cuda.manual_seed_all(seed)","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:32:25.944505Z","iopub.execute_input":"2024-10-06T12:32:25.945346Z","iopub.status.idle":"2024-10-06T12:32:25.950755Z","shell.execute_reply.started":"2024-10-06T12:32:25.945303Z","shell.execute_reply":"2024-10-06T12:32:25.949553Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"set_random_seed(42)\ndevice = torch.device(\"cuda\" if torch.cuda.is_available() else \"cpu\")\nmodel = ST1CNN().to(device)","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:32:26.330728Z","iopub.execute_input":"2024-10-06T12:32:26.331127Z","iopub.status.idle":"2024-10-06T12:32:26.352190Z","shell.execute_reply.started":"2024-10-06T12:32:26.331088Z","shell.execute_reply":"2024-10-06T12:32:26.351157Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train_dataset = DataGetter(df_labels,df_transforms,\"st1\")\ntrain_dataset.load_series_paths(df_split_st1[df_split_st1[\"split\"] == \"train\"][\"study_id\"])\n\nval_dataset = DataGetter(df_labels,df_transforms,\"st1\")\nval_dataset.load_series_paths(df_split_st1[df_split_st1[\"split\"] == \"val\"][\"study_id\"])\n\n# Calculating class weights\nstudies = [int(os.path.normpath(x).split(os.path.sep)[-2]) for x in train_dataset.series_paths]\ntrain_studies = [int(df_labels[df_labels[\"study_id\"] == x][\"colabel\"].item()) for x in studies]\nclass_sample_counts = np.bincount(train_studies)\nclass_weights = 1. / class_sample_counts\nsamples_weights = class_weights[train_studies]\n\n# Creating the sampler for the DataLoader\nsampler = WeightedRandomSampler(samples_weights, len(samples_weights), replacement=False)\n\nbatch_size = 12\nepoch_nb = 30\ntrain_loader = DataLoader(train_dataset, batch_size=batch_size, sampler=sampler, num_workers=0)\nval_loader = DataLoader(val_dataset, batch_size=batch_size, shuffle=True, num_workers=0)\n\noptimizer = optim.Adam(model.parameters(), lr=0.00006)\ncriterions = []\nfor i in range(0,100,4):\n    criterions.append(nn.CrossEntropyLoss(weight=criterion_weights[i:i+4].to(device)))\n\n#optimizer = optim.RMSprop(model.parameters(),lr=0.005, weight_decay=1e-5, centered=True)\ntrain_losses = []\nval_losses = []\nfor iepoch in range(0,epoch_nb):\n    print(f\"--------------Epoch: {iepoch}--------------\")\n    train_loss = 0\n    nb_batch = 0\n    model.train()\n    for serie, labels in train_loader:\n        nb_batch += 1\n        \"\"\"sums = [0,0,0,0]\n        for i in range(0,100,4):\n            sums[0] += torch.sum(labels[:,0,i]).item()\n            sums[1] += torch.sum(labels[:,0,i+1]).item()\n            sums[2] += torch.sum(labels[:,0,i+2]).item()\n            sums[3] += torch.sum(labels[:,0,i+3]).item()\n        print([x/batch_size for x in sums])\"\"\"\n        serie, labels = serie.to(device), labels.squeeze(1).to(device)\n        optimizer.zero_grad()  \n        outputs = model(serie)\n        targets_loss = criterions[0](outputs[:,0:4], labels[:,0:4])\n        for i in range(4,100,4):\n            targets_loss += criterions[i//4](outputs[:,i:i+4], labels[:,i:i+4])\n        targets_loss.backward()\n        optimizer.step()\n        train_loss += targets_loss.item()\n    train_loss = train_loss/(nb_batch*25)\n    train_losses.append(train_loss)\n    print(\"train loss:\",train_loss)\n    model.eval()  \n    with torch.no_grad(): \n        val_loss = 0\n        nb_batch = 0\n        for serie, labels in val_loader:\n            nb_batch += 1\n            serie, labels = serie.to(device), labels.squeeze(1).to(device)\n            outputs = model(serie)\n            targets_loss = criterions[0](outputs[:,0:4], labels[:,0:4])\n            for i in range(4,100,4):\n                targets_loss += criterions[i//4](outputs[:,i:i+4], labels[:,i:i+4])\n            val_loss += targets_loss.item()\n        val_loss = val_loss/(nb_batch*25)\n        val_losses.append(val_loss)\n        print(\"val loss:\",val_loss)","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:33:07.676759Z","iopub.execute_input":"2024-10-06T12:33:07.677171Z","iopub.status.idle":"2024-10-06T12:35:53.239430Z","shell.execute_reply.started":"2024-10-06T12:33:07.677131Z","shell.execute_reply":"2024-10-06T12:35:53.238135Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"torch.save(model, \"../working/model.pth\")","metadata":{"execution":{"iopub.status.busy":"2024-10-06T12:31:49.546253Z","iopub.status.idle":"2024-10-06T12:31:49.546639Z","shell.execute_reply.started":"2024-10-06T12:31:49.546456Z","shell.execute_reply":"2024-10-06T12:31:49.546476Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}