{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "25f736d8",
   "metadata": {},
   "source": [
    "# Using Subject Space ROIs from Freesurfer\n",
    "\n",
    "An example using the AFQ API to find bundles\n",
    "as defined by endpoint ROIs from freesurfer.\n",
    "This example can be modified to work with ROIs\n",
    "in subject space from pipelines other than freesurfer."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7a8a79f6",
   "metadata": {},
   "outputs": [],
   "source": [
    "import os.path as op\n",
    "\n",
    "import nibabel as nib\n",
    "import plotly\n",
    "import numpy as np\n",
    "\n",
    "from AFQ.api.group import GroupAFQ\n",
    "import AFQ.data.fetch as afd\n",
    "from AFQ.definitions.image import RoiImage\n",
    "import AFQ.api.bundle_dict as abd"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b34763d7",
   "metadata": {},
   "source": [
    "## Get some example data\n",
    "\n",
    "Retrieves High angular resolution diffusion imaging (HARDI) dataset from\n",
    "Stanford's Vista Lab\n",
    "\n",
    "See https://purl.stanford.edu/ng782rw8378 for details on dataset.\n",
    "\n",
    "The data for the first subject and first session are downloaded locally\n",
    "(by default into the users home directory) under:\n",
    "\n",
    "`.dipy/stanford_hardi/`\n",
    "\n",
    "Anatomical data (`anat`) and Diffusion-weighted imaging data (`dwi`) are\n",
    "then extracted, formatted to be BIDS compliant, and placed in the AFQ\n",
    "data directory (by default in the users home directory) under:\n",
    "\n",
    "`AFQ_data/stanford_hardi/`\n",
    "\n",
    "This data represents the required preprocessed diffusion data necessary for\n",
    "initializing the GroupAFQ object (which we will do next)\n",
    "\n",
    "The clear_previous_afq is used to remove any previous runs of the afq object\n",
    "stored in the AFQ_data/stanford_hardi/ BIDS directory. Set it to None if\n",
    "you want to use the results of previous runs. Setting it to \"track\"\n",
    "as here will only clear derivatives that depend on the tractography stage\n",
    "(i.e., bundle delination and tract profile calculation),\n",
    "as well as the tractography itself, to save time on recomputation.\n",
    "If you want to only clear derivatives that depend on bundle delineation,\n",
    "and keep the tractography, you can set clear_previous_afq to\n",
    "\"recog\" instead."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "899e17fb",
   "metadata": {},
   "outputs": [],
   "source": [
    "afd.organize_stanford_data(clear_previous_afq=\"track\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4573c7eb",
   "metadata": {},
   "source": [
    "## Generate left thalamus ROI from freesurfer segmentation file\n",
    "1. Load the segmentation file that was generated by Freesurfer for\n",
    "   the specific subject.\n",
    "2. Identify the left thalamus within the file, which has the label\n",
    "   number 41\n",
    "3. Create a Nifti image representing the left thalamus ROI:\n",
    "   - Assign a value of 1 to the voxels that Freesurfer\n",
    "     has labeled as 41 (i.e., the left thalamus).\n",
    "   - Assign a value of 0 to all other voxels.\n",
    "This binary mask format is the expected input for pyAFQ when\n",
    "dealing with subject space ROIs. If it's already in binary format,\n",
    "there is no need to do this step."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "0879f0fd",
   "metadata": {},
   "outputs": [],
   "source": [
    "freesurfer_subject_folder = op.join(\n",
    "    afd.afq_home, \"stanford_hardi\",\n",
    "    \"derivatives\", \"freesurfer\",\n",
    "    \"sub-01\", \"ses-01\",\n",
    "    \"anat\")\n",
    "\n",
    "seg_file = nib.load(op.join(\n",
    "    freesurfer_subject_folder, \"sub-01_ses-01_seg.nii.gz\"))\n",
    "left_thal = seg_file.get_fdata() == 41\n",
    "nib.save(\n",
    "    nib.Nifti1Image(\n",
    "        left_thal.astype(np.float32),\n",
    "        seg_file.affine),\n",
    "    op.join(\n",
    "        freesurfer_subject_folder,\n",
    "        \"sub-01_ses-01_desc-leftThal_mask.nii.gz\"))\n",
    "\n",
    "# Fetch LV1 ROI\n",
    "# which was already generated using the process above\n",
    "afd.fetch_stanford_hardi_lv1()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f47d8138",
   "metadata": {},
   "source": [
    "## Set tractography parameters (optional)\n",
    "We make this tracking_params which we will pass to the GroupAFQ object\n",
    "which specifies that we want 10,000 seeds randomly distributed\n",
    "only within the endpoint ROIs and not throughout the white matter.\n",
    "This is controlled by passing\n",
    "`\"seed_mask\": RoiImage()` in the `tracking_params` dict.\n",
    "\n",
    "We only do this to make this example faster and consume less space."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f40631b1",
   "metadata": {},
   "outputs": [],
   "source": [
    "tracking_params = dict(n_seeds=10000,\n",
    "                       random_seeds=True,\n",
    "                       rng_seed=42,\n",
    "                       seed_mask=RoiImage(use_endpoints=True))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "799f784e",
   "metadata": {},
   "source": [
    "## Define custom `BundleDict` object\n",
    "In a typical `BundleDict` object, ROIs are passed as paths to Nifti files.\n",
    "Here, we define ROIs as dictionaries instead, containing BIDS filters.\n",
    "Then pyAFQ can find the respective ROI for each subject and session."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "693cff24",
   "metadata": {},
   "outputs": [],
   "source": [
    "bundles = abd.BundleDict({\n",
    "    \"L_OR\": {\n",
    "        \"start\": {\n",
    "            \"scope\": \"freesurfer\",\n",
    "            \"suffix\": \"mask\",\n",
    "            \"desc\": \"leftThal\"},\n",
    "        \"end\": {\n",
    "            \"scope\": \"freesurfer\",\n",
    "            \"suffix\": \"anat\",\n",
    "            \"desc\": \"LV1\"\n",
    "        },\n",
    "        \"cross_midline\": False,\n",
    "        \"space\": \"subject\"\n",
    "    }})"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b7ec2ce1",
   "metadata": {},
   "source": [
    "## Initialize a GroupAFQ object:\n",
    "\n",
    "Creates a GroupAFQ object, that encapsulates tractometry,\n",
    "passing in our custom bundle info. Then we run the pipeline\n",
    "and generate a visualization of the bundle we found."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6280680e",
   "metadata": {},
   "outputs": [],
   "source": [
    "myafq = GroupAFQ(\n",
    "    bids_path=op.join(afd.afq_home, 'stanford_hardi'),\n",
    "    dwi_preproc_pipeline='vistasoft',\n",
    "    t1_preproc_pipeline='freesurfer',\n",
    "    tracking_params=tracking_params,\n",
    "    bundle_info=bundles)\n",
    "\n",
    "bundle_html = myafq.export(\"indiv_bundles_figures\")\n",
    "plotly.io.show(bundle_html[\"01\"][\"L_OR\"])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "920c5847",
   "metadata": {},
   "source": [
    ":::{only} html\n",
    "{download}`Download as Jupyter Notebook <use_subject_space_rois_from_freesurfer.ipynb>`\n",
    ":::"
   ]
  }
 ],
 "metadata": {
  "file_format": "mystnb",
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "name": "python",
   "pygments_lexer": "ipython3"
  },
  "mystnb": {
   "execution_mode": "off"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}