{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# An example workflow for voxel-wise encoding models using a BIDS app\n",
    "\n",
    "This shows how to (for a BIDS compliant dataset) extract features, save them in BIDS format, and run a BIDS app for voxel-wise encoding models.\n",
    "We are going to use [this](https://openneuro.org/datasets/ds002322/versions/1.0.4) dataset.\n",
    "\n",
    "#### *Warning*: Executing this notebook will download the full dataset."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!aws s3 sync --no-sign-request s3://openneuro.org/ds002322 ds002322-download/"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Extracting a stimulus representation\n",
    "\n",
    "The dataset in question consists of fMRI activity recorded of several participants while they listened to a reading of the first chapter of Lewis Carroll’s Alice in Wonderland.\n",
    "First we want to extract a stimulus representation that we can use - I chose a Mel spectrogram for demonstration.\n",
    "[This](https://github.com/mjboos/audio2bidsstim/) small Python script extracts such a representation and saves it in a BIDS compliant format.\n",
    "\n",
    "If you get an error that `sndfile library` was not found, you will need to use conda to install it."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import json\n",
    "# these are the parameters for extracting a Mel spectrogram\n",
    "# for computational ease in this example we want 1 sec segments of 31 Mel frequencies with a max frequency of * KHz\n",
    "mel_params = {'n_mels': 31, 'sr': 16000, 'hop_length': 16000, 'n_fft': 16000, 'fmax': 8000}\n",
    "with open('config.json', 'w+') as fl:\n",
    "    json.dump(mel_params, fl)\n",
    "\n",
    "!git clone https://github.com/mjboos/audio2bidsstim/\n",
    "!pip install -r audio2bidsstim/requirements.txt\n",
    "!python audio2bidsstim/wav_files_to_bids_tsv.py ds002322-download/stimuli/DownTheRabbitHoleFinal_mono_exp120_NR16_pad.wav -c config.json\n",
    "!ls -l"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we must copy these files into the BIDS dataset directory according to [these](https://bids-specification.readthedocs.io/en/stable/04-modality-specific-files/06-physiological-and-other-continuous-recordings.html) specifications.\n",
    "We are going to use the `derivatives` folder for the already preprocessed data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!cp DownTheRabbitHoleFinal_mono_exp120_NR16_pad.tsv.gz ds002322-download/derivatives/task-alice_stim.tsv.gz\n",
    "!cp DownTheRabbitHoleFinal_mono_exp120_NR16_pad.json ds002322-download/derivatives/sub-18/sub-18_task-alice_stim.json"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And, lastly, because for this dataset the derivatives folder is missing timing information for the BOLD files - we are only interested in the TR - we have to copy that as well."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!cp ds002322-download/sub-18/sub-18_task-alice_bold.json ds002322-download/derivatives/sub-18/sub-18_task-alice_bold.json "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Running the analysis\n",
    "\n",
    "Now we're all set and can run our encoding analysis. This analysis uses standard Ridge regression, and we're going to specify some additional parameters here."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ridge_params = {'alphas': [1e-1, 1, 100, 1000], 'n_splits': 3, 'normalize': True}\n",
    "\n",
    "# and for lagging the stimulus as well - we want to include 6 sec stimulus segments to predict fMRI\n",
    "lagging_params = {'lag_time': 6}\n",
    "with open('encoding_config.json', 'w+') as fl:\n",
    "    json.dump(ridge_params, fl)\n",
    "    \n",
    "with open('lagging_config.json', 'w+') as fl:\n",
    "    json.dump(lagging_params, fl)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "Now we just need [this](https://github.com/mjboos/voxelwiseencoding) BIDS app for running the analysis.\n",
    "Running this cell will fit voxel-wise encoding models, which right now need about 8 Gig of RAM. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Using Docker to run the voxelwise-encoding BIDS app\n",
    "\n",
    "You can use Docker to build/get an image that already includes all libraries:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!git clone https://github.com/mjboos/voxelwiseencoding\n",
    "!mkdir output\n",
    "# we need to mount a config folder for our json files\n",
    "!mkdir config\n",
    "!cp *config.json config/\n",
    "!docker run -i --rm -v ds002322-download/derivatives:bids_dataset/:ro -v config/:/config:ro -v output/:/output mjboos/voxelwiseencoding /bids_dataset /output --task alice --skip_bids_validator --participant_label 18 --preprocessing-config /config/lagging_config.json --encoding-config /config/encoding_config.json --detrend --standardize zscore "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Alternative: run the module directly\n",
    "\n",
    "Alternatively you can install the required libraries directly and run the Python script yourself."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!git clone https://github.com/mjboos/voxelwiseencoding\n",
    "!pip install -r voxelwiseencoding/requirements.txt\n",
    "!mkdir output\n",
    "!python voxelwiseencoding/run.py ds002322-download/derivatives output --task alice --skip_bids_validator --participant_label 18 --preprocessing-config lagging_config.json --encoding-config encoding_config.json --detrend --standardize zscore "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we'll have some ridge regressions saved in output, as well as scores saved as a Nifti file - which we can visualize.\n",
    "First we load the scores - we have one volume containing the scores per fold - and average them and then plot them via Nilearn."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from nilearn.image import mean_img\n",
    "mean_scores = mean_img('output/sub-18_task-alice_scores.nii.gz')\n",
    "\n",
    "from nilearn import plotting\n",
    "plotting.plot_stat_map(mean_scores, threshold=0.1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And voilà, we see that we can predict activity in the auditory areas."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python (mne)",
   "language": "python",
   "name": "mne"
  },
  "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.7.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}