{"nbformat": 4, "nbformat_minor": 1, "cells": [{"cell_type": "markdown", "source": ["The notebook is an extension of the inofficial IO functions. Here we organize the data a bit better for a standard classification / regression problem. The example creates two DataFrames. One for the input (patient ID and all associated files, and header tags) the other for the output (a single vector sized 17 with each position corresponding to one of the 'zones' mentioned in the overview. The value is then between 0 and 1 corresponding with the probability of that zone. This can thus be directly used as an input to a neural network or random forest for training."], "metadata": {"_cell_guid": "86c3aef9-19ab-4720-b5b8-5b58b292a3ee", "_uuid": "a49a92d4732cae422c49c40cd5fc269dda1397f4", "_execution_state": "idle"}}, {"cell_type": "markdown", "source": ["## Read header"], "metadata": {"_cell_guid": "e69d0f57-8a32-4bfd-b505-2a29e6404f7d", "_uuid": "87d8411dc21ac9edf0fc0f0dbc0f037fd3c9f34a", "_execution_state": "idle"}}, {"cell_type": "code", "outputs": [], "execution_count": null, "source": ["import numpy as np\n", "import os\n", "import matplotlib\n", "import matplotlib.pyplot as plt\n", "from glob import glob\n", "import pandas as pd\n", "import seaborn as sns\n", "plt.style.use('ggplot')\n", "sns.set_style(\"whitegrid\", {'axes.grid' : False})\n", "import matplotlib.animation as mpl_animation\n", "matplotlib.rc('animation', html='html5')\n", "%matplotlib inline"], "metadata": {"collapsed": true, "_cell_guid": "b0bc98f4-2ffb-4e63-bf0d-7b629ec11e25", "_uuid": "6b7df4027deff38111a496a65687f345234204b9", "_execution_state": "idle"}}, {"cell_type": "code", "outputs": [], "execution_count": null, "source": ["def read_header(infile):\n", "    \"\"\"Read image header (first 512 bytes)\n", "    \"\"\"\n", "    h = dict()\n", "    with open(infile, 'r+b') as fid:\n", "        h['filename'] = b''.join(np.fromfile(fid, dtype = 'S1', count = 20))\n", "        h['parent_filename'] = b''.join(np.fromfile(fid, dtype = 'S1', count = 20))\n", "        h['comments1'] = b''.join(np.fromfile(fid, dtype = 'S1', count = 80))\n", "        h['comments2'] = b''.join(np.fromfile(fid, dtype = 'S1', count = 80))\n", "        h['energy_type'] = np.fromfile(fid, dtype = np.int16, count = 1)\n", "        h['config_type'] = np.fromfile(fid, dtype = np.int16, count = 1)\n", "        h['file_type'] = np.fromfile(fid, dtype = np.int16, count = 1)\n", "        h['trans_type'] = np.fromfile(fid, dtype = np.int16, count = 1)\n", "        h['scan_type'] = np.fromfile(fid, dtype = np.int16, count = 1)\n", "        h['data_type'] = np.fromfile(fid, dtype = np.int16, count = 1)\n", "        h['date_modified'] = b''.join(np.fromfile(fid, dtype = 'S1', count = 16))\n", "        h['frequency'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['mat_velocity'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['num_pts'] = np.fromfile(fid, dtype = np.int32, count = 1)\n", "        h['num_polarization_channels'] = np.fromfile(fid, dtype = np.int16, count = 1)\n", "        h['spare00'] = np.fromfile(fid, dtype = np.int16, count = 1)\n", "        h['adc_min_voltage'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['adc_max_voltage'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['band_width'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['spare01'] = np.fromfile(fid, dtype = np.int16, count = 5)\n", "        h['polarization_type'] = np.fromfile(fid, dtype = np.int16, count = 4)\n", "        h['record_header_size'] = np.fromfile(fid, dtype = np.int16, count = 1)\n", "        h['word_type'] = np.fromfile(fid, dtype = np.int16, count = 1)\n", "        h['word_precision'] = np.fromfile(fid, dtype = np.int16, count = 1)\n", "        h['min_data_value'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['max_data_value'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['avg_data_value'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['data_scale_factor'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['data_units'] = np.fromfile(fid, dtype = np.int16, count = 1)\n", "        h['surf_removal'] = np.fromfile(fid, dtype = np.uint16, count = 1)\n", "        h['edge_weighting'] = np.fromfile(fid, dtype = np.uint16, count = 1)\n", "        h['x_units'] = np.fromfile(fid, dtype = np.uint16, count = 1)\n", "        h['y_units'] = np.fromfile(fid, dtype = np.uint16, count = 1)\n", "        h['z_units'] = np.fromfile(fid, dtype = np.uint16, count = 1)\n", "        h['t_units'] = np.fromfile(fid, dtype = np.uint16, count = 1)\n", "        h['spare02'] = np.fromfile(fid, dtype = np.int16, count = 1)\n", "        h['x_return_speed'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['y_return_speed'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['z_return_speed'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['scan_orientation'] = np.fromfile(fid, dtype = np.int16, count = 1)\n", "        h['scan_direction'] = np.fromfile(fid, dtype = np.int16, count = 1)\n", "        h['data_storage_order'] = np.fromfile(fid, dtype = np.int16, count = 1)\n", "        h['scanner_type'] = np.fromfile(fid, dtype = np.int16, count = 1)\n", "        h['x_inc'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['y_inc'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['z_inc'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['t_inc'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['num_x_pts'] = np.fromfile(fid, dtype = np.int32, count = 1)\n", "        h['num_y_pts'] = np.fromfile(fid, dtype = np.int32, count = 1)\n", "        h['num_z_pts'] = np.fromfile(fid, dtype = np.int32, count = 1)\n", "        h['num_t_pts'] = np.fromfile(fid, dtype = np.int32, count = 1)\n", "        h['x_speed'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['y_speed'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['z_speed'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['x_acc'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['y_acc'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['z_acc'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['x_motor_res'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['y_motor_res'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['z_motor_res'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['x_encoder_res'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['y_encoder_res'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['z_encoder_res'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['date_processed'] = b''.join(np.fromfile(fid, dtype = 'S1', count = 8))\n", "        h['time_processed'] = b''.join(np.fromfile(fid, dtype = 'S1', count = 8))\n", "        h['depth_recon'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['x_max_travel'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['y_max_travel'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['elevation_offset_angle'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['roll_offset_angle'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['z_max_travel'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['azimuth_offset_angle'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['adc_type'] = np.fromfile(fid, dtype = np.int16, count = 1)\n", "        h['spare06'] = np.fromfile(fid, dtype = np.int16, count = 1)\n", "        h['scanner_radius'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['x_offset'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['y_offset'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['z_offset'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['t_delay'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['range_gate_start'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['range_gate_end'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['ahis_software_version'] = np.fromfile(fid, dtype = np.float32, count = 1)\n", "        h['spare_end'] = np.fromfile(fid, dtype = np.float32, count = 10)\n", "    return h"], "metadata": {"collapsed": true, "_cell_guid": "b669e9d4-389e-471b-a26e-7b9f59724af4", "_uuid": "d2bbd58d7b33d594459aa6bbb8fc4147668a1ee1", "_execution_state": "idle"}}, {"cell_type": "markdown", "source": ["## Read image data"], "metadata": {"_cell_guid": "32371ee0-c472-4247-8ab7-bb5b12c164ce", "_uuid": "81bc02741351c9f640d41ab18b60663721e82246", "_execution_state": "idle"}}, {"cell_type": "code", "outputs": [], "execution_count": null, "source": ["from collections import namedtuple\n", "from warnings import warn\n", "ScanData = namedtuple('ScanData', ['header', 'data', 'real', 'imag', 'extension'])\n", "def read_data(infile):\n", "    \"\"\"Read any of the 4 types of image files, returns a numpy array of the image contents\n", "    \"\"\"\n", "    _, extension = os.path.splitext(infile)\n", "    sd_dict = {'header': None, 'data': None, 'real': None, 'imag': None, 'extension': extension}\n", "    \n", "    h = read_header(infile)\n", "    sd_dict['header'] = h\n", "    nx = int(h['num_x_pts'])\n", "    ny = int(h['num_y_pts'])\n", "    nt = int(h['num_t_pts'])\n", "    with open(infile, 'rb') as fid:\n", "        fid.seek(512) #skip header\n", "        if extension == '.aps' or extension == '.a3daps':\n", "            if(h['word_type']==7): #float32\n", "                data = np.fromfile(fid, dtype = np.float32, count = nx * ny * nt)\n", "            elif(h['word_type']==4): #uint16\n", "                data = np.fromfile(fid, dtype = np.uint16, count = nx * ny * nt)\n", "            data = data * h['data_scale_factor'] #scaling factor\n", "            data = data.reshape(nx, ny, nt, order='F').copy() #make N-d image\n", "        elif extension == '.a3d':\n", "            if(h['word_type']==7): #float32\n", "                data = np.fromfile(fid, dtype = np.float32, count = nx * ny * nt)\n", "            elif(h['word_type']==4): #uint16\n", "                data = np.fromfile(fid, dtype = np.uint16, count = nx * ny * nt)\n", "            data = data * h['data_scale_factor'] #scaling factor\n", "            data = data.reshape(nx, nt, ny, order='F').copy() #make N-d image\n", "        elif extension == '.ahi':\n", "            data = np.fromfile(fid, dtype = np.float32, count = 2* nx * ny * nt)\n", "            data = data.reshape(2, ny, nx, nt, order='F').copy()\n", "            real = data[0,:,:,:].copy()\n", "            imag = data[1,:,:,:].copy()\n", "            sd_dict['real'] = real\n", "            sd_dict['imag'] = imag\n", "        else:\n", "            warn('Extension not really supported: {}'.format(extension), RuntimeWarning)\n", "            data = None\n", "        sd_dict['data'] = data\n", "        \n", "    return ScanData(**sd_dict)"], "metadata": {"_cell_guid": "220f1039-77ee-4d0a-a7ce-a340f23f961b", "_uuid": "c27f70b9e82bfb7de2329bffb79965586f126d22", "_execution_state": "idle"}}, {"cell_type": "markdown", "source": ["# Plotting Functions\n", "There are a couple of different types of data so we have a couple of different plotting functions\n", "- .aps is for animations\n", "- .a3d is for 3d images"], "metadata": {"_cell_guid": "3542ab26-0712-4caf-9e6b-d88c9e27d925", "_uuid": "d9856629c0f9df190227e19f6d2352571b24d060", "_execution_state": "idle"}}, {"cell_type": "code", "outputs": [], "execution_count": null, "source": ["from skimage.util.montage import montage2d\n", "def plot_montage(sdata):\n", "    if sdata.data is not None:\n", "        print('input data shape', sdata.data.shape)\n", "        fig = plt.figure(figsize = (16,16))\n", "        ax = fig.add_subplot(111)\n", "        ax.imshow(montage2d(np.flipud(sdata.data).swapaxes(0,2)), cmap = 'viridis')\n", "        return fig\n", "def plot_mip(sdata, mip_func = np.max):\n", "    if sdata.data is not None:\n", "        print('input data shape', sdata.data.shape)\n", "        fig, m_axs = plt.subplots(1, 3, figsize = (16,16))\n", "        n_data = np.flipud(sdata.data).swapaxes(0,2)\n", "        for i, (c_name, c_ax)  in enumerate(zip(['xy', 'yz', 'xz'], m_axs)):\n", "            c_ax.imshow(mip_func(n_data,i), cmap = 'viridis')\n", "            c_ax.set_title('%s MIP Projection' % (c_name))\n", "        return fig\n", "            \n", "def plot_animation(sdata):\n", "    if sdata.data is not None:\n", "        print('input data shape', sdata.data.shape)\n", "        fig = plt.figure(figsize = (16,16))\n", "        ax = fig.add_subplot(111)\n", "        def animate(i):\n", "            im = ax.imshow(np.flipud(sdata.data[:,:,i].transpose()), cmap = 'viridis')\n", "            return [im]\n", "        return mpl_animation.FuncAnimation(fig, animate, frames=range(0,sdata.data.shape[2]), interval=200, blit=True)"], "metadata": {"_cell_guid": "8adccb4c-fbdb-4285-a00d-50b9dd710990", "_uuid": "c09380d28fe20fe9d6b0a2af82cea9f9323b28df", "_execution_state": "idle"}}, {"cell_type": "code", "outputs": [], "execution_count": null, "source": ["from keras.utils.np_utils import to_categorical\n", "base_dir = os.path.join('..', 'input')\n", "label_df = pd.read_csv(os.path.join(base_dir, 'stage1_labels.csv'))\n", "label_df['ImageId'] = label_df['Id'].map(lambda x: x.split('_')[0])\n", "label_df['ImageZoneId'] = label_df['Id'].map(lambda x: int(x.split('_')[1][4:]))\n", "# create a vector with each category being an image zone\n", "label_df['ImageZoneVec'] = label_df.apply(\n", "    lambda c_row: [c_row['Probability']*to_categorical(c_row['ImageZoneId']-1, \n", "                                                      num_classes = label_df['ImageZoneId'].max())],1)\n", "label_df\n", "label_df.sample(3)"], "metadata": {"_cell_guid": "4f079efa-d2eb-4909-8065-75268a887533", "_uuid": "548e5734b5892a743ab56ed557dae6353abf63af", "_execution_state": "idle"}}, {"cell_type": "markdown", "source": ["# Create a Aggregate Labelset"], "metadata": {"_cell_guid": "92c73417-9e78-4fcd-a224-f84b9c2f7005", "_uuid": "884c4c2a2206251e6f18fced0ff5105ecba0c21b"}}, {"cell_type": "code", "outputs": [], "execution_count": null, "source": ["def vec_agg(x):\n", "    rslt = dict()\n", "    for col in x.columns:\n", "        rslt[col]=x[col].tolist()\n", "    return pd.Series(rslt)\n", "\n", "agg_label_df = label_df[['ImageId','ImageZoneVec']].groupby('ImageId').apply(vec_agg).drop('ImageId', axis = 1)\n", "agg_label_df = agg_label_df.reset_index()\n", "agg_label_df['ImageZoneVec'] = agg_label_df['ImageZoneVec'].map(lambda x: np.sum(np.hstack(x)[0],0))\n", "agg_label_df['AverageProbability'] = agg_label_df['ImageZoneVec'].map(lambda x: np.mean(x))\n", "agg_label_df.sample(2)"], "metadata": {"_cell_guid": "65f4a507-dfc7-4451-ac12-50f73cc664a3", "_uuid": "f882a7fa9eb54e671fc5f89ea6182ea9afc3942a", "_execution_state": "idle"}}, {"cell_type": "code", "outputs": [], "execution_count": null, "source": ["files_df = pd.DataFrame([dict(ImageId = os.path.splitext(os.path.basename(x))[0], \n", "                              path = x,\n", "                             extension = os.path.splitext(x)[1]) \n", "                         for x in glob(os.path.join(base_dir, 'sample', '*'))])\n", "files_df.sample(2)"], "metadata": {"_cell_guid": "5549b30c-7146-414a-a61c-c62e9e896211", "_uuid": "8a2bd296a5c9e30c09255e08a62a3a1c8068534d", "_execution_state": "idle"}}, {"cell_type": "code", "outputs": [], "execution_count": null, "source": ["def path_agg(x):\n", "    rslt = dict()\n", "    for c_ext, c_path in zip(x['extension'], x['path']):\n", "        rslt[c_ext.replace('.', '')] = c_path\n", "        if c_ext == '.a3d':\n", "            # read the full header for the 3d files\n", "            for i,j in read_header(c_path).items():\n", "                rslt[i] = j\n", "    return pd.Series(rslt)\n", "full_files_df = files_df.groupby('ImageId').apply(path_agg).reset_index()\n", "full_files_df.sample(2)"], "metadata": {"_cell_guid": "a630a60b-e9b4-46f4-88da-f2306af17f98", "_uuid": "dc9004eec093f906c0984fe895d5d0ca7126ed43", "_execution_state": "idle"}}, {"cell_type": "code", "outputs": [], "execution_count": null, "source": ["if False:\n", "    for _, c_row in full_files_df.sample(1).iterrows():\n", "        full_3d = read_data(c_row['a3d'])\n", "        fig = plot_mip(full_3d, np.sum)"], "metadata": {"scrolled": false, "_cell_guid": "1808e004-7389-4517-aa4d-af656b40c282", "_uuid": "077a6f6bc5264ac68b4550dbeaa180fa2620f6b2", "_execution_state": "idle"}}, {"cell_type": "code", "outputs": [], "execution_count": null, "source": ["if False:\n", "    for _, c_row in full_files_df.sample(1).iterrows():\n", "        full_aps = read_data(c_row['aps'])\n", "        plot_montage(full_aps)"], "metadata": {"_cell_guid": "fbbea62a-2c66-43c1-9ffb-82dde4976667", "_uuid": "b898efcb49d82aaaa7563dfce43a77a4d5ad5558", "_execution_state": "idle"}}, {"cell_type": "markdown", "source": ["# Combine the labels to the images\n", "Probably makes more sense after analysis has been done"], "metadata": {"_cell_guid": "1685d66c-e102-418c-9ef0-1ce6045c5047", "_uuid": "7317d24d802cdbb6a43531d07500874b47648d57"}}, {"cell_type": "code", "outputs": [], "execution_count": null, "source": ["comb_df = pd.merge(full_files_df, agg_label_df)\n", "comb_df.sample(2)"], "metadata": {"_cell_guid": "77b7207b-2cb9-4ea2-a8f6-ae4b7ef6f969", "_uuid": "088f05cbdea06e722de2a6bd65fa118b4439bb89", "_execution_state": "idle"}}, {"cell_type": "code", "outputs": [], "execution_count": null, "source": [], "metadata": {"_cell_guid": "1ea09d30-585b-464d-a09a-448a2e07d0e0", "_uuid": "04e37bb3b7bdebf97178c9d57a8b633311d54d9a", "_execution_state": "idle"}}, {"cell_type": "code", "outputs": [], "execution_count": null, "source": [], "metadata": {"_cell_guid": "ca09da93-94d2-475f-889f-b117bafba01d", "_uuid": "0945342a5de6fb60b1cc0763c813c2b7f7a5242c", "_execution_state": "idle"}}, {"cell_type": "code", "outputs": [], "execution_count": null, "source": [], "metadata": {"collapsed": true, "_cell_guid": "35df424a-4ccf-49f9-8c04-0f30421a6ece", "_uuid": "5daebc86182c658ae83ed882bc591bc69b2ba9b4", "_execution_state": "idle"}}, {"cell_type": "code", "outputs": [], "execution_count": null, "source": [], "metadata": {"_cell_guid": "267fa396-88b4-4d63-98f7-90d7ca8a79d6", "_uuid": "8295e0591aa3913e1381761ba0555e7bb0a19fe1", "_execution_state": "idle"}}, {"cell_type": "code", "outputs": [], "execution_count": null, "source": [], "metadata": {"collapsed": true, "_cell_guid": "cf0f5256-1891-4979-a6cf-ea43043ac5f8", "_uuid": "e6aaa56db72089f38de9ce0a17b9e73ec776590d", "_execution_state": "idle"}}], "metadata": {"anaconda-cloud": {}, "language_info": {"nbconvert_exporter": "python", "name": "python", "version": "3.6.3", "mimetype": "text/x-python", "codemirror_mode": {"name": "ipython", "version": 3}, "file_extension": ".py", "pygments_lexer": "ipython3"}, "kernelspec": {"display_name": "Python 3", "language": "python", "name": "python3"}}}