pipeline_2a_functionalsessions_preprocessing.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Generate BIDS-formatted folder structure for fmriprep\n",
    "\n",
    "Note that we use a modified version of the T1w-scans as input to ensure skull stripping within fmriprep works.\n",
    "\n",
    "\n",
    "Pipeline:\n",
    "1. Create BIDS-formatted folder structure for fmriprep, ensuring that there is a single T1w-scan per subject (ses-me, otherwise ses-se)\n",
    "2. Multiply T1w image by $\\frac{inv2}{inv2+\\gamma}$. Gamma is currently set to 100 (determined experimentally). This procedure removes the background noise while minimizing the amount of remaining bias in the image. See https://www.ncbi.nlm.nih.gov/pmc/articles/PMC4059664/ for rationale behind the used equation.\n",
    "4. Run fmriprep"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os\n",
    "import numpy as np\n",
    "import shutil\n",
    "import glob\n",
    "import nibabel as nib\n",
    "from nilearn import plotting, masking\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 62,
   "metadata": {},
   "outputs": [],
   "source": [
    "raw_dir = './data/raw/bids'\n",
    "fmriprep_raw_dir = './data/raw/bids_fmriprep'\n",
    "os.makedirs(fmriprep_raw_dir, exist_ok=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 63,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "sub-01, ses-me... ses-se... \n",
      "sub-02, ses-me... ses-se... \n",
      "sub-03, ses-me... ses-se... \n",
      "sub-04, ses-me... ses-se... \n",
      "sub-05, ses-me... ses-se... \n",
      "sub-06, ses-me... ses-se... \n",
      "sub-07, ses-me... ses-se... \n",
      "sub-08, ses-me... ses-se... \n",
      "sub-09, ses-me... ses-se... \n",
      "sub-10, ses-me... ses-se... \n",
      "sub-11, ses-me... ses-se... \n",
      "sub-12, ses-me... ses-se... \n",
      "sub-13, ses-me... ses-se... \n",
      "sub-14, ses-me... ses-se... \n",
      "sub-15, ses-me... ses-se... \n",
      "sub-16, ses-me... ses-se... \n",
      "sub-17, ses-me... ses-se... \n",
      "sub-18, ses-me... ses-se... "
     ]
    }
   ],
   "source": [
    "fns = os.listdir(raw_dir)\n",
    "fns.sort()\n",
    "fns = [os.path.join(raw_dir, fn) for fn in fns]\n",
    "\n",
    "for fn in fns:\n",
    "    # files can be just copied\n",
    "    if os.path.isfile(fn):\n",
    "        this_file = fn.split('/')[-1]\n",
    "        shutil.copy2(fn, os.path.join(fmriprep_raw_dir, this_file))\n",
    "    else:\n",
    "        subj_id = fn.split('/')[-1]\n",
    "        print()\n",
    "        print(subj_id, end=', ')\n",
    "        \n",
    "        subj_dir = os.path.join(fmriprep_raw_dir, subj_id)\n",
    "        os.makedirs(subj_dir, exist_ok=True)\n",
    "        \n",
    "        for ses in ['ses-me', 'ses-se']:\n",
    "            print(ses, end='... ')\n",
    "            # copy func, fmap\n",
    "            if os.path.isdir(os.path.join(fn, ses)):\n",
    "                shutil.copytree(os.path.join(fn, ses), os.path.join(subj_dir, ses))\n",
    "        \n",
    "        if os.path.isdir(os.path.join(fn, 'ses-me', 'anat')):\n",
    "            # ses-me has anatomy, remove ses-ses anatomy if this also exists\n",
    "            if os.path.isdir(os.path.join(fn, 'ses-se', 'anat')):\n",
    "                shutil.rmtree(os.path.join(subj_dir, 'ses-se', 'anat'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 64,
   "metadata": {},
   "outputs": [],
   "source": [
    "def t1w_reduce_noise(t1w_path, inv2_path, gamma=100):\n",
    "    ### OVERWRITES T1w images!\n",
    "    t1w = nib.load(t1w_path)\n",
    "    inv2 = nib.load(inv2_path)\n",
    "    \n",
    "    t1w_data = t1w.get_data()\n",
    "    inv2_data = inv2.get_data()\n",
    "    \n",
    "    t1w_data = t1w_data * (inv2_data / (inv2_data+gamma))\n",
    "    new_t1w = nib.Nifti1Image(t1w_data, t1w.affine)\n",
    "    anat_dir = '/'.join(t1w_path.split('/')[:-1])\n",
    "    os.system('touch {}/.t1w_modified'.format(anat_dir))\n",
    "    nib.save(new_t1w, t1w_path)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 65,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "01, 02, 03, 04, 05, 06, 07, 08, 09, 10, 11, 12, 13, 14, 15, 16, 17, 18, "
     ]
    }
   ],
   "source": [
    "# generate \"new\" T1w-images\n",
    "fns = os.listdir(fmriprep_raw_dir)\n",
    "fns.sort()\n",
    "fns = [os.path.join(fmriprep_raw_dir, fn) for fn in fns]\n",
    "cmds = []\n",
    "\n",
    "for fn in fns:\n",
    "    if os.path.isfile(fn):\n",
    "        continue\n",
    "    else:\n",
    "        subj_id = fn.split('/')[-1].split('-')[-1]\n",
    "        print(subj_id, end=', ')\n",
    "        \n",
    "        if os.path.isdir(os.path.join(fn, 'ses-me', 'anat')):\n",
    "            # ses-me has anatomy\n",
    "            anat_location = os.path.join(fn, 'ses-me', 'anat')\n",
    "            anat_ses = 'me'\n",
    "        else:\n",
    "            # anatomy is from ses-se\n",
    "            anat_location = os.path.join(fn, 'ses-se', 'anat')\n",
    "            anat_ses = 'se'\n",
    "\n",
    "        # generate \"T1w\" img\n",
    "        if os.path.exists(os.path.join(anat_location, '.t1w_modified')):\n",
    "            continue\n",
    "        t1w_path = os.path.join(anat_location, 'sub-{}_ses-{}_T1w.nii'.format(subj_id, anat_ses))\n",
    "        inv2_path = os.path.join(anat_location, 'sub-{}_ses-{}_inv-2_MPRAGE.nii'.format(subj_id, anat_ses))\n",
    "        t1w_reduce_noise(t1w_path, inv2_path, gamma=100)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "T1w_imgs = glob.glob('./data/raw/bids_fmriprep/sub-*/ses-*/anat/*_T1w.nii.gz')\n",
    "T1w_imgs.sort()\n",
    "\n",
    "for T1 in T1w_imgs:\n",
    "    print(T1)\n",
    "    plotting.plot_anat(T1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Afterwards, we run fmriprep. Note that I used a slightly modified version of fmriprep 1.2.6. That is, since fmriprep 1.2, it automatically generates an optimal combination of the multiple echos in an ME dataset. I did not want this, rather, I wanted it to simply process each echo separately (except for HMC - which works better with the first echo only).\n",
    "\n",
    "Fmriprep was run on the LISA system in Amsterdam. The docker image used is `stevenm/mefmri_1.2.6.2`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "/usr/local/miniconda/lib/python3.7/site-packages/requests/__init__.py:91: RequestsDependencyWarning: urllib3 (1.24.1) or chardet (3.0.4) doesn't match a supported version!\n",
      "  RequestsDependencyWarning)\n",
      "/usr/local/miniconda/lib/python3.7/site-packages/requests/__init__.py:91: RequestsDependencyWarning: urllib3 (1.24.1) or chardet (3.0.4) doesn't match a supported version!\n",
      "  RequestsDependencyWarning)\n",
      "/usr/local/miniconda/lib/python3.7/importlib/_bootstrap.py:219: ImportWarning: can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "  return f(*args, **kwds)\n",
      "/usr/local/miniconda/lib/python3.7/importlib/_bootstrap.py:219: ImportWarning: can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "  return f(*args, **kwds)\n",
      "/usr/local/miniconda/lib/python3.7/importlib/_bootstrap.py:219: ImportWarning: can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "  return f(*args, **kwds)\n",
      "/usr/local/miniconda/lib/python3.7/site-packages/nilearn/datasets/neurovault.py:16: DeprecationWarning: Using or importing the ABCs from 'collections' instead of from 'collections.abc' is deprecated, and in 3.8 it will stop working\n",
      "  from collections import Container\n",
      "/usr/local/miniconda/lib/python3.7/site-packages/skimage/__init__.py:80: ResourceWarning: unclosed file <_io.TextIOWrapper name='/usr/local/miniconda/lib/python3.7/site-packages/pytest.py' mode='r' encoding='utf-8'>\n",
      "  imp.find_module('pytest')\n",
      "/usr/local/miniconda/lib/python3.7/importlib/_bootstrap.py:219: ImportWarning: can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "  return f(*args, **kwds)\n",
      "/usr/local/miniconda/lib/python3.7/site-packages/networkx/classes/reportviews.py:95: DeprecationWarning: Using or importing the ABCs from 'collections' instead of from 'collections.abc' is deprecated, and in 3.8 it will stop working\n",
      "  from collections import Mapping, Set, Iterable\n",
      "/usr/local/miniconda/lib/python3.7/site-packages/bids/layout/bids_layout.py:121: ResourceWarning: unclosed file <_io.TextIOWrapper name='/home/neuro/work_dir/data/raw/bids_fmriprep/dataset_description.json' mode='r' encoding='UTF-8'>\n",
      "  self.description = json.load(open(target, 'r'))\n",
      "/usr/local/miniconda/lib/python3.7/site-packages/grabbit/core.py:436: ResourceWarning: unclosed file <_io.TextIOWrapper name='/usr/local/miniconda/lib/python3.7/site-packages/bids/layout/config/bids.json' mode='r' encoding='UTF-8'>\n",
      "  domain = json.load(open(domain, 'r'))\n",
      "/usr/local/miniconda/lib/python3.7/site-packages/requests/__init__.py:91: RequestsDependencyWarning: urllib3 (1.24.1) or chardet (3.0.4) doesn't match a supported version!\n",
      "  RequestsDependencyWarning)\n",
      "/usr/local/miniconda/lib/python3.7/importlib/_bootstrap.py:219: ImportWarning: can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "  return f(*args, **kwds)\n",
      "/usr/local/miniconda/lib/python3.7/importlib/_bootstrap.py:219: ImportWarning: can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "  return f(*args, **kwds)\n",
      "/usr/local/miniconda/lib/python3.7/importlib/_bootstrap.py:219: ImportWarning: can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "  return f(*args, **kwds)\n",
      "/usr/local/miniconda/lib/python3.7/importlib/_bootstrap.py:219: ImportWarning: can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "  return f(*args, **kwds)\n",
      "/usr/local/miniconda/lib/python3.7/site-packages/nilearn/datasets/neurovault.py:16: DeprecationWarning: Using or importing the ABCs from 'collections' instead of from 'collections.abc' is deprecated, and in 3.8 it will stop working\n",
      "  from collections import Container\n",
      "/usr/local/miniconda/lib/python3.7/site-packages/skimage/__init__.py:80: ResourceWarning: unclosed file <_io.TextIOWrapper name='/usr/local/miniconda/lib/python3.7/site-packages/pytest.py' mode='r' encoding='utf-8'>\n",
      "  imp.find_module('pytest')\n",
      "/usr/local/miniconda/lib/python3.7/importlib/_bootstrap.py:219: ImportWarning: can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "  return f(*args, **kwds)\n",
      "[WARNING] This document format requires a nonempty <title> element.\n",
      "  Please specify either 'title' or 'pagetitle' in the metadata.\n",
      "  Falling back to 'CITATION'\n",
      "Captured warning (<class 'ImportWarning'>): can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "Captured warning (<class 'ImportWarning'>): can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "Captured warning (<class 'ImportWarning'>): can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "Captured warning (<class 'ImportWarning'>): can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "Captured warning (<class 'DeprecationWarning'>): Using or importing the ABCs from 'collections' instead of from 'collections.abc' is deprecated, and in 3.8 it will stop working\n",
      "Captured warning (<class 'ResourceWarning'>): unclosed file <_io.TextIOWrapper name='/usr/local/miniconda/lib/python3.7/site-packages/pytest.py' mode='r' encoding='utf-8'>\n",
      "Captured warning (<class 'ImportWarning'>): can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "Captured warning (<class 'PendingDeprecationWarning'>): the matrix subclass is not the recommended way to represent matrices or deal with linear algebra (see https://docs.scipy.org/doc/numpy/user/numpy-for-matlab-users.html). Please adjust your code to use regular ndarray.\n",
      "/usr/local/miniconda/lib/python3.7/site-packages/requests/__init__.py:91: RequestsDependencyWarning: urllib3 (1.24.1) or chardet (3.0.4) doesn't match a supported version!\n",
      "  RequestsDependencyWarning)\n",
      "/usr/local/miniconda/lib/python3.7/importlib/_bootstrap.py:219: ImportWarning: can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "  return f(*args, **kwds)\n",
      "/usr/local/miniconda/lib/python3.7/importlib/_bootstrap.py:219: ImportWarning: can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "  return f(*args, **kwds)\n",
      "/usr/local/miniconda/lib/python3.7/importlib/_bootstrap.py:219: ImportWarning: can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "  return f(*args, **kwds)\n",
      "Captured warning (<class 'PendingDeprecationWarning'>): the matrix subclass is not the recommended way to represent matrices or deal with linear algebra (see https://docs.scipy.org/doc/numpy/user/numpy-for-matlab-users.html). Please adjust your code to use regular ndarray.\n",
      "/usr/local/miniconda/lib/python3.7/site-packages/nilearn/datasets/neurovault.py:16: DeprecationWarning: Using or importing the ABCs from 'collections' instead of from 'collections.abc' is deprecated, and in 3.8 it will stop working\n",
      "  from collections import Container\n",
      "/usr/local/miniconda/lib/python3.7/site-packages/requests/__init__.py:91: RequestsDependencyWarning: urllib3 (1.24.1) or chardet (3.0.4) doesn't match a supported version!\n",
      "  RequestsDependencyWarning)\n",
      "/usr/local/miniconda/lib/python3.7/importlib/_bootstrap.py:219: ImportWarning: can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "  return f(*args, **kwds)\n",
      "/usr/local/miniconda/lib/python3.7/importlib/_bootstrap.py:219: ImportWarning: can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "  return f(*args, **kwds)\n",
      "/usr/local/miniconda/lib/python3.7/importlib/_bootstrap.py:219: ImportWarning: can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "  return f(*args, **kwds)\n",
      "/usr/local/miniconda/lib/python3.7/site-packages/nilearn/datasets/neurovault.py:16: DeprecationWarning: Using or importing the ABCs from 'collections' instead of from 'collections.abc' is deprecated, and in 3.8 it will stop working\n",
      "  from collections import Container\n",
      "/usr/local/miniconda/lib/python3.7/site-packages/requests/__init__.py:91: RequestsDependencyWarning: urllib3 (1.24.1) or chardet (3.0.4) doesn't match a supported version!\n",
      "  RequestsDependencyWarning)\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "/usr/local/miniconda/lib/python3.7/site-packages/requests/__init__.py:91: RequestsDependencyWarning: urllib3 (1.24.1) or chardet (3.0.4) doesn't match a supported version!\n",
      "  RequestsDependencyWarning)\n",
      "/usr/local/miniconda/lib/python3.7/importlib/_bootstrap.py:219: ImportWarning: can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "  return f(*args, **kwds)\n",
      "/usr/local/miniconda/lib/python3.7/importlib/_bootstrap.py:219: ImportWarning: can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "  return f(*args, **kwds)\n",
      "/usr/local/miniconda/lib/python3.7/importlib/_bootstrap.py:219: ImportWarning: can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "  return f(*args, **kwds)\n",
      "/usr/local/miniconda/lib/python3.7/site-packages/nilearn/datasets/neurovault.py:16: DeprecationWarning: Using or importing the ABCs from 'collections' instead of from 'collections.abc' is deprecated, and in 3.8 it will stop working\n",
      "  from collections import Container\n",
      "/usr/local/miniconda/lib/python3.7/site-packages/requests/__init__.py:91: RequestsDependencyWarning: urllib3 (1.24.1) or chardet (3.0.4) doesn't match a supported version!\n",
      "  RequestsDependencyWarning)\n",
      "/usr/local/miniconda/lib/python3.7/site-packages/requests/__init__.py:91: RequestsDependencyWarning: urllib3 (1.24.1) or chardet (3.0.4) doesn't match a supported version!\n",
      "  RequestsDependencyWarning)\n",
      "/usr/local/miniconda/lib/python3.7/importlib/_bootstrap.py:219: ImportWarning: can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "  return f(*args, **kwds)\n",
      "/usr/local/miniconda/lib/python3.7/importlib/_bootstrap.py:219: ImportWarning: can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "  return f(*args, **kwds)\n",
      "/usr/local/miniconda/lib/python3.7/importlib/_bootstrap.py:219: ImportWarning: can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "  return f(*args, **kwds)\n",
      "/usr/local/miniconda/lib/python3.7/site-packages/nilearn/datasets/neurovault.py:16: DeprecationWarning: Using or importing the ABCs from 'collections' instead of from 'collections.abc' is deprecated, and in 3.8 it will stop working\n",
      "  from collections import Container\n",
      "/usr/local/miniconda/lib/python3.7/site-packages/skimage/__init__.py:80: ResourceWarning: unclosed file <_io.TextIOWrapper name='/usr/local/miniconda/lib/python3.7/site-packages/pytest.py' mode='r' encoding='utf-8'>\n",
      "  imp.find_module('pytest')\n",
      "/usr/local/miniconda/lib/python3.7/importlib/_bootstrap.py:219: ImportWarning: can't resolve package from __spec__ or __package__, falling back on __name__ and __path__\n",
      "  return f(*args, **kwds)\n"
     ]
    }
   ],
   "source": [
    "# !fmriprep ./data/raw/bids_fmriprep ./data/deriv/fmriprep_126_run2 participant \\\n",
    "# --participant-label 01 \\\n",
    "# -w ./workflow_folders \\\n",
    "# --n_cpus 10 \\\n",
    "# --mem-mb 100000 \\\n",
    "# --anat-only \\\n",
    "# -v \\\n",
    "# --output-space T1w template \\\n",
    "# --template-resampling-grid native > $(date +\"run_%Y%m%d_%H%M%S.txt\")"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}