{"cells":[{"metadata":{"_kg_hide-input":true,"_kg_hide-output":true,"trusted":true},"cell_type":"code","source":"!pip install pynrrd","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"from tqdm import tqdm\nimport os\nfrom random import randint\n\nimport numpy as np\nimport pandas as pd\n\nimport nibabel as nib\nimport pydicom as pdm\nimport nilearn as nl\nimport nilearn.plotting as nlplt\nimport nrrd\nimport h5py\n\nimport matplotlib.pyplot as plt\nfrom matplotlib import cm\nimport matplotlib.animation as anim\n\nimport imageio\nfrom skimage.transform import resize\nfrom skimage.util import montage\n\nfrom IPython.display import Image as show_gif\n\nimport warnings\nwarnings.simplefilter(\"ignore\")","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"class ImageToGIF:\n    \"\"\"Create GIF without saving image files.\"\"\"\n    def __init__(self,\n                 size=(600, 400), \n                 xy_text=(80, 10),\n                 dpi=100, \n                 cmap='CMRmap'):\n\n        self.fig = plt.figure()\n        self.fig.set_size_inches(size[0] / dpi, size[1] / dpi)\n        self.xy_text = xy_text\n        self.cmap = cmap\n        \n        self.ax = self.fig.add_axes([0, 0, 1, 1])\n        self.ax.set_xticks([])\n        self.ax.set_yticks([])\n        self.images = []\n \n    def add(self, *args, label, with_mask=True):\n        \n        image = args[0]\n        mask = args[-1]\n        plt.set_cmap(self.cmap)\n        plt_img = self.ax.imshow(image, animated=True)\n        if with_mask:\n            plt_mask = self.ax.imshow(np.ma.masked_where(mask == False, mask),\n                                      alpha=0.7, animated=True)\n\n        plt_text = self.ax.text(*self.xy_text, label, color='red')\n        to_plot = [plt_img, plt_mask, plt_text] if with_mask else [plt_img, plt_text]\n        self.images.append(to_plot)\n        plt.close()\n \n    def save(self, filename, fps):\n        animation = anim.ArtistAnimation(self.fig, self.images)\n        animation.save(filename, writer='imagemagick', fps=fps)\n        \n        \nclass Image3dToGIF3d:\n    \"\"\"\n    Displaying 3D images in 3d axes.\n    Parameters:\n        img_dim: shape of cube for resizing.\n        figsize: figure size for plotting in inches.\n    \"\"\"\n    def __init__(self, \n                 img_dim: tuple = (55, 55, 55),\n                 figsize: tuple = (15, 10),\n                ):\n        \"\"\"Initialization.\"\"\"\n        self.img_dim = img_dim\n        print(img_dim)\n        self.figsize = figsize\n    \n    def _explode(self, data: np.ndarray):\n        \"\"\"\n        Takes: array and return an array twice as large in each dimension,\n        with an extra space between each voxel.\n        \"\"\"\n        shape_arr = np.array(data.shape)\n        size = shape_arr[:3] * 2 - 1\n        exploded = np.zeros(np.concatenate([size, shape_arr[3:]]),\n                            dtype=data.dtype)\n        exploded[::2, ::2, ::2] = data\n        return exploded\n\n    def _expand_coordinates(self, indices: np.ndarray):\n        x, y, z = indices\n        x[1::2, :, :] += 1\n        y[:, 1::2, :] += 1\n        z[:, :, 1::2] += 1\n        return x, y, z\n    \n    def _normalize(self, arr: np.ndarray):\n        \"\"\"Normilize image value between 0 and 1.\"\"\"\n        return arr / arr.max()\n    \n    def _scale_by(self, arr: np.ndarray, factor: int):\n        \"\"\"\n        Scale 3d Image to factor.\n        Parameters:\n            arr: 3d image for scalling.\n            factor: factor for scalling.\n        \"\"\"\n        mean = np.mean(arr)\n        return (arr - mean) * factor + mean\n    \n    def get_transformed_data(self, data: np.ndarray):\n        \"\"\"Data transformation: normalization, scaling, resizing.\"\"\"\n        norm_data = np.clip(self._normalize(data)-0.1, 0, 1) ** 0.4\n        scaled_data = np.clip(self._scale_by(norm_data, 2) - 0.1, 0, 1)\n        resized_data = resize(scaled_data, self.img_dim, mode='constant')\n        return resized_data\n    \n    def plot_cube(self,\n                  cube,\n                  title: str = '', \n                  init_angle: int = 0,\n                  make_gif: bool = False,\n                  path_to_save: str = 'filename.gif'\n                 ):\n        \"\"\"\n        Plot 3d data.\n        Parameters:\n            cube: 3d data\n            title: title for figure.\n            init_angle: angle for image plot (from 0-360).\n            make_gif: if True create gif from every 5th frames from 3d image plot.\n            path_to_save: path to save GIF file.\n            \"\"\"\n        cube = self._normalize(cube)\n\n        facecolors = cm.gist_stern(cube)\n        facecolors[:,:,:,-1] = cube\n        facecolors = self._explode(facecolors)\n\n        filled = facecolors[:,:,:,-1] != 0\n        x, y, z = self._expand_coordinates(np.indices(np.array(filled.shape) + 1))\n\n        with plt.style.context(\"dark_background\"):\n\n            fig = plt.figure(figsize=self.figsize)\n            ax = fig.gca(projection='3d')\n\n            ax.view_init(30, init_angle)\n            ax.set_xlim(right = self.img_dim[0] * 2)\n            ax.set_ylim(top = self.img_dim[1] * 2)\n            ax.set_zlim(top = self.img_dim[2] * 2)\n            ax.set_title(title, fontsize=18, y=1.05)\n\n            ax.voxels(x, y, z, filled, facecolors=facecolors, shade=False)\n\n            if make_gif:\n                images = []\n                for angle in tqdm(range(0, 360, 5)):\n                    ax.view_init(30, angle)\n                    fname = str(angle) + '.png'\n\n                    plt.savefig(fname, dpi=120, format='png', bbox_inches='tight')\n                    images.append(imageio.imread(fname))\n                    #os.remove(fname)\n                imageio.mimsave(path_to_save, images)\n                plt.close()\n\n            else:\n                plt.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Nifti File"},{"metadata":{"trusted":true},"cell_type":"code","source":"sample_filename = '../input/brats20-dataset-training-validation/BraTS2020_TrainingData/MICCAI_BraTS2020_TrainingData/BraTS20_Training_001/BraTS20_Training_001_flair.nii'\nsample_filename_mask = '../input/brats20-dataset-training-validation/BraTS2020_TrainingData/MICCAI_BraTS2020_TrainingData/BraTS20_Training_001/BraTS20_Training_001_seg.nii'\n\nsample_img = nib.load(sample_filename)\nsample_img = np.asanyarray(sample_img.dataobj)\nsample_mask = nib.load(sample_filename_mask)\nsample_mask = np.asanyarray(sample_mask.dataobj)\nprint(\"img shape ->\", sample_img.shape)\nprint(\"mask shape ->\", sample_mask.shape)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"visualization of 3d data in 2d slices"},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"slice_n = 100\nfig, ax = plt.subplots(2, 3, figsize=(25, 15))\n\nax[0, 0].imshow(sample_img[slice_n, :, :])\nax[0, 0].set_title(f\"image slice number {slice_n} along the x-axis\", fontsize=18, color=\"red\")\nax[1, 0].imshow(sample_mask[slice_n, :, :])\nax[1, 0].set_title(f\"mask slice {slice_n} along the x-axis\", fontsize=18, color=\"red\")\n\nax[0, 1].imshow(sample_img[:, slice_n, :])\nax[0, 1].set_title(f\"image slice number {slice_n} along the y-axis\", fontsize=18, color=\"red\")\nax[1, 1].imshow(sample_mask[:, slice_n, :])\nax[1, 1].set_title(f\"mask slice number {slice_n} along the y-axis\", fontsize=18, color=\"red\")\n\nax[0, 2].imshow(sample_img[:, :, slice_n])\nax[0, 2].set_title(f\"image slice number {slice_n} along the z-axis\", fontsize=18, color=\"red\")\nax[1, 2].imshow(sample_mask[:, :, slice_n])\nax[1, 2].set_title(f\"mask slice number {slice_n}along the z-axis\", fontsize=18, color=\"red\")\nfig.tight_layout()\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"## matching colormaps\n#Greys_r RdGy_r  CMRmap afmhot binary_r bone copper cubehelix gist_heat gist_stern gnuplot hot inferno magma nipy_spectral\n\nsample_data_gif = ImageToGIF()\nlabel = sample_filename.replace('/', '.').split('.')[-2]\nfilename = f'{label}_3d_2d.gif'\n\nfor i in range(sample_img.shape[0]):\n    image = np.rot90(sample_img[i])\n    mask = np.clip(np.rot90(sample_mask[i]), 0, 1)\n    sample_data_gif.add(image, mask, label=f'{label}_{str(i)}')\n \nsample_data_gif.save(filename, fps=15)\nshow_gif(filename, format='png')","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"image = np.rot90(montage(sample_img))\nmask = np.rot90(montage(sample_mask)) \nmask = np.clip(mask, 0, 1)\n\nfig, ax1 = plt.subplots(1, 1, figsize = (20, 20))\nax1.imshow(image, cmap ='bone')\nax1.imshow(np.ma.masked_where(mask == False, mask),\n           cmap='cool', alpha=0.6, animated=True)\nfig.savefig(f'{label}_3d_to_2d.png', format='png', bbox_inches='tight', pad_iches=0.0)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"and 3d"},{"metadata":{"trusted":true},"cell_type":"code","source":"%%time\ntitle = sample_filename.replace(\".\", \"/\").split(\"/\")[-2]\nfilename = title+\"_3d.gif\"\n\ndata_to_3dgif = Image3dToGIF3d()#img_dim = (120, 120, 78)\ntransformed_data = data_to_3dgif.get_transformed_data(sample_img)\ndata_to_3dgif.plot_cube(\n    transformed_data[:38, :47, :35],#[:77, :105, :55]\n    title=title,\n    make_gif=True,\n    path_to_save=filename\n)\nshow_gif(filename, format='png')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"visualization with nilearn"},{"metadata":{"trusted":true},"cell_type":"code","source":"niimg = nl.image.load_img(sample_filename)\nnimask = nl.image.load_img(sample_filename_mask)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"fig, axes = plt.subplots(nrows=4, figsize=(30, 40))\n\n\nnlplt.plot_anat(niimg,\n                title='BraTS20_Training_001_flair.nii plot_anat',\n                axes=axes[0])\n\nnlplt.plot_epi(niimg,\n               title='BraTS20_Training_001_flair.nii plot_epi',\n               axes=axes[1])\n\nnlplt.plot_img(niimg,\n               title='BraTS20_Training_001_flair.nii plot_img',\n               axes=axes[2])\n\nnlplt.plot_roi(nimask, \n               title='BraTS20_Training_001_flair.nii with mask plot_roi',\n               bg_img=niimg, \n               axes=axes[3], cmap='Paired')\n\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"t1_niimg  = nl.image.load_img('../input/brats20-dataset-training-validation/BraTS2020_TrainingData/MICCAI_BraTS2020_TrainingData/BraTS20_Training_001/BraTS20_Training_001_t1.nii')\nt2_niimg  = nl.image.load_img('../input/brats20-dataset-training-validation/BraTS2020_TrainingData/MICCAI_BraTS2020_TrainingData/BraTS20_Training_001/BraTS20_Training_001_t2.nii')\n\nfig, axes = plt.subplots(nrows=2, figsize=(30, 20))\n\nnlplt.plot_epi(t1_niimg, title=\"BraTS20_Training_001_t1.nii plot_epi\", axes=axes[0])\nnlplt.plot_epi(t2_niimg, title=\"BraTS20_Training_001_t2.nii plot_epi\", axes=axes[1])\nplt.show()","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# DICOM File"},{"metadata":{"trusted":true},"cell_type":"code","source":"sample_path = '../input/osic-pulmonary-fibrosis-progression/train/ID00015637202177877247924'\nsample_path_files = sorted(os.listdir(sample_path), key=lambda x: int(x[:-4]))","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"visualization of 3d data in 2d slices"},{"metadata":{"trusted":true},"cell_type":"code","source":"sample_data_gif = ImageToGIF(size=(768, 768),\n                             xy_text=(250, 15))\n\nlabel = sample_path.split('/')[-1]\nfor i in range(len(sample_path_files)):\n    path = os.path.join(sample_path, sample_path_files[i])\n    image = pdm.dcmread(path).pixel_array\n    sample_data_gif.add(image, label=f'{label}_{str(i)}', with_mask=False)\n \nsample_data_gif.save(f'{label}.gif', fps=15)\nshow_gif(f'{label}.gif', format='png')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# Nrrd File"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"def read_nrrd_file(path: str, \n                   tensor_shape: tuple ) -> np.ndarray:\n    if os.path.exists(path):\n        tensor = nrrd.read(path)[0]                             \n        tensor = np.flip(tensor, -1)                   # Warning! slice order of images and masks does not match.\n    else: \n        tensor = np.zeros(tensor_shape, dtype=np.float32)\n    return tensor\n\n\ndef nrrd_to_numpy(id_: str, tensor_shape: tuple):\n    '''\n    Returns:  all id masks in single numpy tensor.\n    '''\n    lung_file_path = '../input/ct-lung-heart-trachea-segmentation/nrrd_lung/nrrd_lung/' + id_ + '_lung.nrrd'\n    heart_file_path  = '../input/ct-lung-heart-trachea-segmentation/nrrd_heart/nrrd_heart/' + id_ + '_heart.nrrd'         # Here path hardcoded\n    trachea_file_path = '../input/ct-lung-heart-trachea-segmentation/nrrd_trachea/nrrd_trachea/' + id_ + '_trachea.nrrd'\n\n    lung_tensor = read_nrrd_file(lung_file_path, tensor_shape)\n    heart_tensor = read_nrrd_file(heart_file_path, tensor_shape)\n    trachea_tensor = read_nrrd_file(trachea_file_path, tensor_shape)\n    \n    # chek if all tensors  have the same shape.\n    if not (lung_tensor.shape == heart_tensor.shape == trachea_tensor.shape):\n        #print(lung_tensor.shape, heart_tensor.shape, trachea_tensor.shape)\n        print(\"problem with id:\", id_)\n        return \n        \n    # now each tensor channel is a mask with a unique label\n    full_mask = np.stack([lung_tensor, heart_tensor, trachea_tensor])\n\n    # reorient the axes from CHWB to BWHC\n    full_mask = np.moveaxis(full_mask,\n                            [0, 1, 2, 3],\n                            [3, 2, 1, 0]).astype(np.float32)\n\n    return full_mask","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"id_ = 'ID00015637202177877247924'\n\nsample_masks = nrrd_to_numpy(id_, (768, 768))\nsample_masks.shape","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"visualization of 3d data in 2d slices"},{"metadata":{"trusted":true},"cell_type":"code","source":"sample_data_gif = ImageToGIF(size=(768, 768),\n                             xy_text=(250, 15))\n\nlabel = sample_path.split('/')[-1] + '_mask'\nfor i in range(sample_masks.shape[0]):\n    sample_data_gif.add(sample_masks[i],label=f'{label}_{str(i)}', with_mask=False)\n \nsample_data_gif.save(f'{label}.gif', fps=15)\nshow_gif(f'{label}.gif', format='png')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"and images with masks"},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"class ImageToGIF:\n    \"\"\"Create GIF without saving image files.\"\"\"\n    def __init__(self,\n                 size=(768, 768), \n                 xy_text=(250, 15),\n                 dpi=100):\n\n        self.fig = plt.figure()\n        self.fig.set_size_inches(size[0] / dpi, size[1] / dpi)\n        self.xy_text = xy_text\n        \n        self.ax = self.fig.add_axes([0, 0, 1, 1])\n        self.ax.set_xticks([])\n        self.ax.set_yticks([])\n        self.images = []\n \n    def add(self, image, mask, label, with_mask=True):\n        \n        lung, heart, trachea = [mask[:, :, i] for i in range(3)]\n        plt_img = self.ax.imshow(image, cmap=\"bone\", animated=True)\n\n        # Overlaying segmentation masks\n        plt_mask1 = self.ax.imshow(np.ma.masked_where(lung == False, lung),\n                               cmap='cool', alpha=0.3, animated=True)\n        plt_mask2 = self.ax.imshow(np.ma.masked_where(heart == False, heart),\n                               cmap='autumn', alpha=0.3, animated=True)\n        plt_mask3 = self.ax.imshow(np.ma.masked_where(trachea == False, trachea),\n                               cmap='autumn_r', alpha=0.3, animated=True) #cool_r - blue\n        plt_text = self.ax.text(*self.xy_text, label, color='red')\n        to_plot = [plt_img, plt_mask1, plt_mask2, plt_mask3, plt_text] \n        self.images.append(to_plot)\n        plt.close()\n \n    def save(self, filename, fps):\n        animation = anim.ArtistAnimation(self.fig, self.images)\n        animation.save(filename, writer='imagemagick', fps=fps)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":false},"cell_type":"code","source":"sample_data_gif = ImageToGIF()\n\nlabel = sample_path.split('/')[-1] + '_with_masks'\nfor i in range(sample_masks.shape[0]):\n    path = os.path.join(sample_path, sample_path_files[i])\n    image = pdm.dcmread(path).pixel_array\n    mask = sample_masks[i]\n    sample_data_gif.add(image, mask, label=f'{label}_{str(i)}',)\n \nsample_data_gif.save(f'{label}.gif', fps=15)\nshow_gif(f'{label}.gif', format='png')","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"# MATLAB File"},{"metadata":{"trusted":true},"cell_type":"code","source":"sample_filename = '../input/trends-assessment-prediction/fMRI_train/10001.mat'\nmatlab_file = h5py.File(sample_filename)\nprint(matlab_file.keys())\nprint(matlab_file.values())\nprint(matlab_file['SM_feature'][()].shape)","execution_count":null,"outputs":[]},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"def get_spatial_maps(path: str) -> np.ndarray:\n    \"\"\"\n    Read the given .mat file with h5py\n    Switch second and fourth axes (53, 52, 63, 53) -> (53, 53, 63, 52)\n    Make identical ZYX axes (53, 53, 63, 52) -> (53, 63, 63, 63)\n    \n    Params:\n        path - of the MATLAB file\n    \n    Returns:\n        mew_spatial_maps - nympy array, shape = (53, 63, 63, 63)\n        masked_spatial_data : numpy array, shape = (53, 58869)\n        numpy array after filling ZYX dimensions.\n    \"\"\"\n    matlab_file = h5py.File(path)\n\n    spatial_maps = np.moveaxis(matlab_file['SM_feature'][()], [0, 1, 2, 3], [3, 2, 1, 0])\n\n\n    max_lengh_ax = np.argmax([spatial_maps.shape[0], spatial_maps.shape[1], spatial_maps.shape[2]])\n    new_spatial_maps = np.full((*[spatial_maps.shape[max_lengh_ax]]*3, spatial_maps.shape[-1],), 0, dtype='<f8')\n\n    new_spatial_maps[:spatial_maps.shape[0],\n                     :spatial_maps.shape[1],\n                     :spatial_maps.shape[2],\n                     :spatial_maps.shape[3]] = spatial_maps\n    \n    return new_spatial_maps","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"sample_spatial_maps_file = '../input/trends-assessment-prediction/fMRI_train/10009.mat'\nsample_mask_file = '../input/trends-assessment-prediction/fMRI_mask.nii'\ngifname = sample_spatial_maps_file.replace('.', '/').split(\"/\")[-2]\n\nfmri_mask = nl.image.load_img(sample_mask_file)\nfmri_mask = np.asarray(fmri_mask.dataobj)\nspatial_maps = get_spatial_maps(sample_spatial_maps_file)\n\nprint(\"spatial_maps ->\", spatial_maps.shape)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"visualization of 3d data in 2d slices"},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"PATH_TO_SAVE = \"test\"\nif not os.path.exists(PATH_TO_SAVE):\n    os.mkdir(PATH_TO_SAVE)\n    print(f'Folder {PATH_TO_SAVE} created.')","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"num_slices = spatial_maps.shape[0]\nstep = randint(0, 52)\n\nfig, ax = plt.subplots(nrows=1, ncols=3, figsize=(20, 12.5))\nplt.set_cmap('bwr')\n[axi.set_axis_off() for axi in ax.ravel()]\n\nfor slice_ in np.arange(num_slices):\n    ax[0].imshow(np.rot90(spatial_maps[slice_, :, :, step], axes=(-2,-1)))\n    ax[1].imshow(np.rot90(spatial_maps[:, slice_, :, step], axes=(-2,-1)))\n    ax[2].imshow(np.rot90(spatial_maps[:, :, slice_, step], axes=(-2,-1)))\n\n    ax[0].set_title('X-axis {}'.format(slice_), fontsize=15,weight='bold', color='white')\n    ax[1].set_title('Y-axis {}'.format(slice_), fontsize=15, weight='bold', color='white')\n    ax[2].set_title('Z-axis {}'.format(slice_), fontsize=15, weight='bold', color='white')\n\n    plt.subplots_adjust(wspace=0, hspace=0)\n    fig.savefig(f'{PATH_TO_SAVE}/{slice_}.png', format='png', bbox_inches='tight', pad_iches=0.0, facecolor='#0e1111') \n    plt.close()","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"images = []\nfilenames = sorted(os.listdir('test'), key=lambda x: int(x[:-4]))\nfor filename in filenames:\n    filename = 'test/' + filename\n    images.append(imageio.imread(filename))\nimageio.mimsave(f'{gifname}.gif', images, duration=0.5)\n\nshow_gif(f'{gifname}.gif', format='png')","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"_kg_hide-output":true,"trusted":true},"cell_type":"code","source":"!rm -r test","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"and 3d "},{"metadata":{"trusted":true,"_kg_hide-input":true},"cell_type":"code","source":"spatial_maps = matlab_file['SM_feature'][()]\nn_map = 22\ntitle = sample_filename.split(\"/\")[-1]\npath_to_save = title + \"_3d.gif\"\n\nX,Y,Z = [], [], [] \nroi = []\n\nfor z, image in enumerate(spatial_maps[n_map,:,:,:]):\n    xx, yy = np.meshgrid(np.linspace(0, image.shape[1],image.shape[1]),\n                         np.linspace(0, image.shape[0],image.shape[0]))\n\n    zz = np.ones(xx.shape) * z\n    xx = xx[[image != 0][0]]\n    yy = yy[[image != 0][0]]\n    zz = np.ones(xx.shape) * z\n    X += list(xx)\n    Y += list(yy)\n    Z += list(zz)\n    roi += list(image[[image != 0][0]])\n\nfig = plt.figure(figsize=(12, 10))\nax = fig.gca(projection='3d')\nax.set_title(title, fontsize=15, y=1.01)\nax.scatter(X, Y, Z, c = plt.cm.bwr(roi), cmap ='bwr', s=1, alpha=1)\n\nimages = []\nfor angle in tqdm(range(0, 360, 5)):\n    ax.view_init(30, angle)\n    fname = str(angle) + '.png'\n    #print(fname)\n    plt.savefig(fname, dpi=120, format='png', bbox_inches='tight')\n    images.append(imageio.imread(fname))\n    #os.remove(fname)\nimageio.mimsave(path_to_save, images)\nplt.close()\n\nshow_gif(path_to_save)","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"visualization using nilearn "},{"metadata":{"trusted":true,"_kg_hide-input":true,"_kg_hide-output":true},"cell_type":"code","source":"!wget https://github.com/Chaogan-Yan/DPABI/raw/master/Templates/ch2better.nii","execution_count":null,"outputs":[]},{"metadata":{"trusted":true},"cell_type":"code","source":"# Loading reference image\nfmri_mask = nl.image.load_img('../input/trends-assessment-prediction/fMRI_mask.nii')\n\n# Reorienting the axis of 3D spatial map\nspatial_maps = np.moveaxis(matlab_file['SM_feature'][()], [0, 1, 2, 3], [3, 2, 1, 0]) \n\n# Loading 3D spatial maps\nspatial_maps_niimg = nl.image.new_img_like(ref_niimg=fmri_mask,\n                                           data=spatial_maps,\n                                           affine=fmri_mask.affine,\n                                           copy_header=True)","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"PATH_TO_SAVE = \"test\"\nif not os.path.exists(PATH_TO_SAVE):\n    os.mkdir(PATH_TO_SAVE)\n    print(f'Folder {PATH_TO_SAVE} created.')","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"for i, img in enumerate(list(nl.image.iter_img(spatial_maps_niimg))):\n    fig, ax = plt.subplots(figsize=(18, 5.5))\n    nlplt.plot_stat_map(stat_map_img=img,\n                bg_img='ch2better.nii',\n                title=f'10001.mat Spatial Map {str(i)} plot_stat_map',\n                axes=ax,\n                threshold=1,\n                display_mode='ortho',\n                annotate=False,\n                draw_cross=True,\n                colorbar=False)\n    plt.subplots_adjust(wspace=0, hspace=0)\n    fig.savefig(f'{PATH_TO_SAVE}/{i}.png', format='png', bbox_inches='tight', pad_iches=0.0, facecolor='#0e1111') \n    plt.close()","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":false,"trusted":true},"cell_type":"code","source":"images = []\nfilenames = sorted(os.listdir('test'), key=lambda x: int(x[:-4]))\nfor filename in filenames:\n    filename = 'test/' + filename\n    images.append(imageio.imread(filename))\nimageio.mimsave(f'{gifname}_nil.gif', images, duration=0.5)\n\nshow_gif(f'{gifname}_nil.gif', format='png')","execution_count":null,"outputs":[]},{"metadata":{"_kg_hide-input":true,"trusted":true},"cell_type":"code","source":"!rm -r test","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"and 3d"},{"metadata":{"trusted":true},"cell_type":"code","source":"img = img = list(nl.image.iter_img(spatial_maps_niimg))[0]\nview = nlplt.view_img_on_surf(img,\n                              title=f'10009.mat Spatial Map 0 view_img_on_surf',\n                              title_fontsize=20,\n                              threshold=1,\n                              black_bg=False)\nview.open_in_browser()\nview","execution_count":null,"outputs":[]},{"metadata":{},"cell_type":"markdown","source":"References:\n+ https://nilearn.github.io/auto_examples/index.html\n+ https://terbium.io/2017/12/matplotlib-3d/\n+ https://tomroelandts.com/articles/how-to-create-animated-gifs-with-python\n+ https://www.kaggle.com/gunesevitan/trends-neuroimaging-data-analysis-3d-features#5.-Component-Spatial-Maps\n+ https://www.kaggle.com/mawanda/fmri-3d-brain-images-dynamic-drawing"}],"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat":4,"nbformat_minor":4}