{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.10.12","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"none","dataSources":[{"sourceId":91249,"databundleVersionId":11294684,"sourceType":"competition"}],"dockerImageVersionId":30918,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"In this notebook, I intend to provide some useful functions for getting started and some helpful visualizations for understanding the data and the task at hand. I hope you find them useful :).","metadata":{}},{"cell_type":"code","source":"import os\nimport cv2\nimport numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true,"execution":{"iopub.status.busy":"2025-03-09T18:14:41.567755Z","iopub.execute_input":"2025-03-09T18:14:41.568407Z","iopub.status.idle":"2025-03-09T18:14:41.574696Z","shell.execute_reply.started":"2025-03-09T18:14:41.568357Z","shell.execute_reply":"2025-03-09T18:14:41.572914Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"#Reading and taking a look at the data\ndata = pd.read_csv('/kaggle/input/byu-locating-bacterial-flagellar-motors-2025/train_labels.csv')\nprint('It looks like the organizers are using -1.0 as a placeholder for when there are no flagellum coordinates in an image.')\ndata.head()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-09T18:14:41.576579Z","iopub.execute_input":"2025-03-09T18:14:41.576972Z","iopub.status.idle":"2025-03-09T18:14:41.622907Z","shell.execute_reply.started":"2025-03-09T18:14:41.576941Z","shell.execute_reply":"2025-03-09T18:14:41.621251Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"#The df names are a little obtuse. Let's change these to be a little more intuitive.\ndata.rename(columns={'Motor axis 0': 'flagellum_z',\n                    'Motor axis 1': 'flagellum_y',\n                    'Motor axis 2': 'flagellum_x',\n                    'Array shape (axis 0)': 'z_dim',\n                    'Array shape (axis 1)': 'y_dim',\n                    'Array shape (axis 2)': 'x_dim'\n                        }, inplace=True)\ndata.head()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-09T18:14:41.625751Z","iopub.execute_input":"2025-03-09T18:14:41.626123Z","iopub.status.idle":"2025-03-09T18:14:41.643563Z","shell.execute_reply.started":"2025-03-09T18:14:41.626094Z","shell.execute_reply":"2025-03-09T18:14:41.642023Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def get_tomogram(directory_path, experiment_name, scaling='None'):\n    \"\"\"\n    Creates a tomogram array (z, x, y) from a path to the directory and name of the experiment.\n\n    Args:\n        directory_path (str): The path to the directory that holds the individual experiments.\n        experiment_name (str): The name of the directory that corresponds to a single experiment, filled with .json files\n        scaling (str): 'None' -> do not scale the tomogram. 'Z' -> standardize the tomogram. 'Max' -> minmax scale the tomogram. Defaults to 'None'.\n\n    Returns:\n        numpy.ndarray: A NumPy array corresponding to the tomogram with shape (z, x, y).\n    \"\"\"\n    path = os.path.join(directory_path, experiment_name)\n    \n    file_list = os.listdir(path)\n    \n    file_list = sorted(file_list, key=lambda x: x)\n    \n    arrays = []\n    for file in file_list:\n        arr = cv2.imread(os.path.join(path, file), cv2.IMREAD_GRAYSCALE)\n\n        arrays.append(arr)\n\n    tomogram = np.stack(arrays, axis=0)\n\n    scaling = str.lower(scaling)\n\n    if scaling == 'z':\n        tomogram = (tomogram-np.mean(tomogram))/(np.maximum(np.std(tomogram),1e-6))\n\n    elif scaling == 'max':\n        tomogram = (tomogram - np.min(tomogram))/np.max(tomogram)\n\n    else:\n        print('Tomogram was not scaled. Available methods for scaling are Z or max.')\n\n    return tomogram","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-09T18:14:41.645932Z","iopub.execute_input":"2025-03-09T18:14:41.646297Z","iopub.status.idle":"2025-03-09T18:14:41.670282Z","shell.execute_reply.started":"2025-03-09T18:14:41.646261Z","shell.execute_reply":"2025-03-09T18:14:41.668542Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def tomogram_histogram(tomogram):\n    \"\"\"\n    Visualize the distribution of pixel intensities in the 3D tomogram.\n\n    Args:\n        tomogram (numpy.ndarray): The (z, x, y) array representing the tomogram data.\n    \"\"\"\n    plt.hist(tomogram.flatten(), color='#AEC6CF', edgecolor='black')\n    plt.xlabel(\"Pixel Intensities\")\n    plt.ylabel(\"Frequency\")\n    plt.title(\"Distribution of Pixel Intensities\")\n    plt.show()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-09T18:14:41.672036Z","iopub.execute_input":"2025-03-09T18:14:41.672555Z","iopub.status.idle":"2025-03-09T18:14:41.692973Z","shell.execute_reply.started":"2025-03-09T18:14:41.672496Z","shell.execute_reply":"2025-03-09T18:14:41.691535Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def view_slice(tomogram, z):\n    \"\"\"\n    Visualize the tomogram at a specified z-axis slice.\n\n    Args:\n        tomogram (numpy.ndarray): The (z, x, y) array representing the tomogram data.\n        z (int): The index of the z-axis slice you would like to visualize.\n    \"\"\"\n\n    image = tomogram[z]\n\n    # Show the image using Matplotlib\n    plt.imshow(image, cmap=\"gray\")\n    plt.axis(\"off\")\n    plt.show()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-09T18:14:41.694451Z","iopub.execute_input":"2025-03-09T18:14:41.694975Z","iopub.status.idle":"2025-03-09T18:14:41.724373Z","shell.execute_reply.started":"2025-03-09T18:14:41.694941Z","shell.execute_reply":"2025-03-09T18:14:41.723079Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def view_slice_circle(tomogram, z, y, x, voxel_spacing, thickness=2, ax=None):\n    \"\"\"\n    Takes a tomogram, coordinates, and voxel spacing, then shows an image where the flagellum is circled.\n\n    Args:\n        tomogram (numpy array): The 3D tomogram.\n        z (float): The Z position of the flagellum to visualize.\n        y (float): The Y position of the flagellum to visualize.\n        x (float): The X position of the flagellum to visualize.\n        voxel_spacing (float): The voxel scaling factor in angstroms.\n        thickness (int, optional): Line thickness of the circle. Defaults to 2. Use -1 to fill the circle.\n        ax (matplotlib.axes.Axes, optional): The axis to plot on. If None, a new figure is created.\n    \"\"\"\n    x, y, z = int(x), int(y), int(z)\n    center = (x, y)\n\n    # Scale the radius based on voxel spacing (nearest angstrom)\n    radius = int(1000 / voxel_spacing)\n\n    if all(v >= 0 for v in [x, y, z]):\n        # Draw the circle on the image\n        image = tomogram[z].copy()\n        color = int(np.max(image))  # Use max intensity for visibility\n        cv2.circle(image, center, radius, color, thickness)\n\n        if ax is None:\n            plt.figure(figsize=(5, 5))\n            plt.imshow(image, cmap=\"gray\")\n            plt.axis(\"off\")\n            plt.show()\n        else:\n            ax.imshow(image, cmap=\"gray\")\n            ax.axis(\"off\")\n    else:\n        if ax:\n            ax.axis(\"off\")  # Hide empty plots in grid\n        print('There is no flagellum in this tomogram. Use \"view_slice\" if you want to see a specific slice of this tomogram.')\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-09T18:14:41.725572Z","iopub.execute_input":"2025-03-09T18:14:41.726000Z","iopub.status.idle":"2025-03-09T18:14:41.749078Z","shell.execute_reply.started":"2025-03-09T18:14:41.725972Z","shell.execute_reply":"2025-03-09T18:14:41.747839Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"def plot_all_flagellum(tomogram, experiment, data):\n    exp_data = data[data['tomo_id'] == experiment].copy()\n    num_images = len(exp_data)\n    cols = 3  # Number of images per row\n    rows = (num_images + cols - 1) // cols  # Compute number of rows\n    \n    fig, axes = plt.subplots(rows, cols, figsize=(cols * 5, rows * 5))\n    axes = np.array(axes).reshape(rows, cols)  # Ensure a 2D array of axes\n    \n    for idx, (i, row) in enumerate(exp_data.iterrows()):\n        z, y, x, spacing = row['flagellum_z'], row['flagellum_y'], row['flagellum_x'], row['Voxel spacing']\n        ax = axes[idx // cols, idx % cols]\n        view_slice_circle(tomogram, z, y, x, spacing, thickness=2, ax=ax)\n        ax.set_title(f\"Flagellum {idx}\")\n\n    # Hide unused subplots\n    for idx in range(num_images, rows * cols):\n        fig.delaxes(axes[idx // cols, idx % cols])\n\n    plt.tight_layout()\n    plt.show()\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-09T18:14:41.750276Z","iopub.execute_input":"2025-03-09T18:14:41.750621Z","iopub.status.idle":"2025-03-09T18:14:41.775051Z","shell.execute_reply.started":"2025-03-09T18:14:41.750580Z","shell.execute_reply":"2025-03-09T18:14:41.773664Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"directory = '/kaggle/input/byu-locating-bacterial-flagellar-motors-2025/train'\nexperiment = 'tomo_00e463' # <-- Change this if you want to look at a different tomogram.\ntom = get_tomogram(directory, experiment, scaling='None')\ntomogram_histogram(tom)\nprint('''This data fits a normal distribution pretty well.\nIt would be good to look through a few more to see if any distributions are impacted by artifacts/outliers.\nConsider experimenting with different scaling techniques to help model training.''')","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-09T18:14:41.778666Z","iopub.execute_input":"2025-03-09T18:14:41.779055Z","iopub.status.idle":"2025-03-09T18:14:53.477644Z","shell.execute_reply.started":"2025-03-09T18:14:41.779027Z","shell.execute_reply":"2025-03-09T18:14:53.476528Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"#Let's take a loook at one of the examples with a flagellum\nplot_all_flagellum(tom, experiment, data)\nprint('''For this competition, it looks like we need to regress the voxel indices where the flagellum meets the body of the cell.\nSurprisingly to me, there can be multiple flagella attached to a single cell. I learned something new today!\nIt also looks like the 1000 angstrom radius we have to predict within is fairly generous.''')","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-09T18:14:53.479098Z","iopub.execute_input":"2025-03-09T18:14:53.479468Z","iopub.status.idle":"2025-03-09T18:14:55.144896Z","shell.execute_reply.started":"2025-03-09T18:14:53.479438Z","shell.execute_reply":"2025-03-09T18:14:55.143799Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"experiments = data[['tomo_id',\n                    'z_dim',\n                    'y_dim',\n                    'x_dim',\n                    'Voxel spacing',\n                    'Number of motors']].groupby('tomo_id').agg('max')","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-09T18:14:55.145760Z","iopub.execute_input":"2025-03-09T18:14:55.146160Z","iopub.status.idle":"2025-03-09T18:14:55.156434Z","shell.execute_reply.started":"2025-03-09T18:14:55.146125Z","shell.execute_reply":"2025-03-09T18:14:55.155301Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"num_flag = experiments['Number of motors'].tolist()\nplt.hist(num_flag, color='#AEC6CF', edgecolor='black')\n# Labels and title\nplt.xlabel(\"Number of Flagellum\")\nplt.ylabel(\"Frequency\")\nplt.title(\"Distribution of Flagellum Per Tomogram\")\nplt.show()\nprint('''A bit more than half the tomograms have a flagellum present.\nthink about this from the perspective of class imbalance.\nOn one hand, our training data is pretty well balanced from a binary classification perspective.\nOn the other hand, the # of pixels that actually contain a flagellum are certainly a minority class compared to background.''')","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-09T18:14:55.157561Z","iopub.execute_input":"2025-03-09T18:14:55.158272Z","iopub.status.idle":"2025-03-09T18:14:55.406421Z","shell.execute_reply.started":"2025-03-09T18:14:55.158233Z","shell.execute_reply":"2025-03-09T18:14:55.405261Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"spacing = experiments['Voxel spacing'].tolist()\nplt.hist(spacing, color='#AEC6CF', edgecolor='black')\n# Labels and title\nplt.xlabel(\"Voxel Spacing\")\nplt.ylabel(\"Frequency\")\nplt.title(\"Voxel spacing per tomogram\")\nplt.show()\nprint('''This plot gives us hints about what distribution of data our models need to be able to handle.\nUnderstanding how to manage this is likely to be important because we do not have this data available at test time.\nRemember, we don't know the test distribution and need to use this to infer what the test distribution is likely to be.''')","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-09T18:14:55.407617Z","iopub.execute_input":"2025-03-09T18:14:55.408025Z","iopub.status.idle":"2025-03-09T18:14:55.784947Z","shell.execute_reply.started":"2025-03-09T18:14:55.407990Z","shell.execute_reply":"2025-03-09T18:14:55.783725Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"experiments['Depth'] = experiments['z_dim']*experiments['Voxel spacing']","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-09T18:14:55.786054Z","iopub.execute_input":"2025-03-09T18:14:55.786474Z","iopub.status.idle":"2025-03-09T18:14:55.792154Z","shell.execute_reply.started":"2025-03-09T18:14:55.786432Z","shell.execute_reply":"2025-03-09T18:14:55.790823Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"num_flag = experiments['Depth'].tolist()\nplt.hist(num_flag, color='#AEC6CF', edgecolor='black')\n# Labels and title\nplt.xlabel(\"Tomogram Depth (Angstroms)\")\nplt.ylabel(\"Frequency\")\nplt.title(\"Distribution of Tomogram Depth\")\nplt.show()\nprint('''It looks like the depth of these experiments varies by a bit less than 4x from min to max and follows a bimodal distribution.\nDoes this variation lead to differeces in pixel intensities? Might be something to look into.\nI leave this to the reader to explore a bit more.''')","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-09T18:14:55.793277Z","iopub.execute_input":"2025-03-09T18:14:55.793663Z","iopub.status.idle":"2025-03-09T18:14:56.020869Z","shell.execute_reply.started":"2025-03-09T18:14:55.793628Z","shell.execute_reply":"2025-03-09T18:14:56.019758Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"experiments['aspect_ratio'] = experiments['x_dim']/experiments['y_dim']\nnum_flag = experiments['aspect_ratio'].tolist()\nplt.hist(num_flag, color='#AEC6CF', edgecolor='black')\n# Labels and title\nplt.xlabel(\"Aspect Ratio\")\nplt.ylabel(\"Frequency\")\nplt.title(\"Distribution of Aspect Ratios (X/Y)\")\nplt.show()\nprint('''For the most part, the X/Y aspect ratios are pretty close to 1.0, but there are a few odd ones.''')","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-09T18:14:56.022175Z","iopub.execute_input":"2025-03-09T18:14:56.022576Z","iopub.status.idle":"2025-03-09T18:14:56.214311Z","shell.execute_reply.started":"2025-03-09T18:14:56.022535Z","shell.execute_reply":"2025-03-09T18:14:56.213163Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"experiments['aspect_ratio'] = experiments['x_dim']/experiments['z_dim']\nnum_flag = experiments['aspect_ratio'].tolist()\nplt.hist(num_flag, color='#AEC6CF', edgecolor='black')\n# Labels and title\nplt.xlabel(\"Aspect Ratio\")\nplt.ylabel(\"Frequency\")\nplt.title(\"Distribution of Aspect Ratios (X/Z)\")\nplt.show()\nprint('''Most of the X/Z ratios are about 3, but there are some others to consider.''')","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-09T18:14:56.215421Z","iopub.execute_input":"2025-03-09T18:14:56.215735Z","iopub.status.idle":"2025-03-09T18:14:56.444959Z","shell.execute_reply.started":"2025-03-09T18:14:56.215669Z","shell.execute_reply":"2025-03-09T18:14:56.443573Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"experiments['aspect_ratio'] = experiments['y_dim']/experiments['z_dim']\nnum_flag = experiments['aspect_ratio'].tolist()\nplt.hist(num_flag, color='#AEC6CF', edgecolor='black')\n# Labels and title\nplt.xlabel(\"Aspect Ratio\")\nplt.ylabel(\"Frequency\")\nplt.title(\"Distribution of Aspect Ratios (Y/Z)\")\nplt.show()\nprint('''Like we saw above, the ratio of Y/Z is typically ~3, but sometimes much lower.''')","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-03-09T18:14:56.446237Z","iopub.execute_input":"2025-03-09T18:14:56.446540Z","iopub.status.idle":"2025-03-09T18:14:56.704747Z","shell.execute_reply.started":"2025-03-09T18:14:56.446513Z","shell.execute_reply":"2025-03-09T18:14:56.703599Z"}},"outputs":[],"execution_count":null}]}