ISC_Analysis.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "69361490",
   "metadata": {},
   "source": [
    "## Load Package"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "076dab50",
   "metadata": {},
   "outputs": [],
   "source": [
    "from isc_standalone import (isc, isfc, bootstrap_isc, permutation_isc,\n",
    "                            timeshift_isc, phaseshift_isc,\n",
    "                            compute_summary_statistic, load_images,\n",
    "                            load_boolean_mask, mask_images,\n",
    "                            MaskedMultiSubjectData)\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import seaborn as sns\n",
    "from scipy.stats import norm, pearsonr, zscore\n",
    "from scipy.spatial.distance import squareform\n",
    "from statsmodels.stats.multitest import multipletests\n",
    "import nibabel as nib\n",
    "\n",
    "from isc_cli import load_data"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4566f69e",
   "metadata": {},
   "source": [
    "## Load Dataset\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "f59c7526",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "data0path = '/Users/gangxinli/Desktop/Internship/Neuro/Neuro_ISC/Data/11Sep/vtc_tal/sub-sid000007_task-movie_run-05_bold_SCCAI_3DMCTS_THPGLMF3c_256_trilin_2x1.nii'\n",
    "data1path = '/Users/gangxinli/Desktop/Internship/Neuro/Neuro_ISC/Data/11Sep/vtc_tal/sub-sid000009_task-movie_run-05_bold_SCCAI_3DMCTS_THPGLMF3c_256_trilin_2x1.nii'\n",
    "data2path = '/Users/gangxinli/Desktop/Internship/Neuro/Neuro_ISC/Data/11Sep/vtc_tal/sub-sid000010_task-movie_run-05_bold_SCCAI_3DMCTS_THPGLMF3c_256_trilin_2x1.nii'\n",
    "data3path = '/Users/gangxinli/Desktop/Internship/Neuro/Neuro_ISC/Data/11Sep/vtc_tal/sub-sid000005_task-movie_run-05_bold_SCCAI_3DMCTS_THPGLMF3c_256_trilin_2x1.nii'\n",
    "mask = '/Users/gangxinli/Desktop/Internship/Neuro/Neuro_ISC/Data/11Sep/vtc_tal/VTC_TALmask.nii'\n",
    "data, affine, header, input_fns = load_data([data0path,data1path,data2path,data3path],mask=mask)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "af833c6c",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Create simple simulated data with high intersubject correlation\n",
    "n_subjects = 4\n",
    "n_TRs = 300\n",
    "n_voxels = 1000\n",
    "def simulated_timeseries(n_subjects, n_TRs, n_voxels=1, noise=1):\n",
    "    signal = np.random.randn(n_TRs, n_voxels // 100)\n",
    "    data = [zscore(np.repeat(signal, 100, axis=1) +\n",
    "                 np.random.randn(n_TRs, n_voxels) * noise,\n",
    "                 axis=0)\n",
    "          for subject in np.arange(n_subjects)]\n",
    "    return data\n",
    "\n",
    "# List of subject datasets\n",
    "data_simulated = simulated_timeseries(n_subjects, n_TRs, n_voxels=n_voxels)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "71923f8c",
   "metadata": {},
   "source": [
    "## Insepect the data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "8c590720",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Sample of the data shape\n",
      "(4, 300, 1000)\n",
      "n_subjects, n_TRs, n_voxels\n"
     ]
    }
   ],
   "source": [
    "print('Sample of the data shape')\n",
    "print(np.array(data_simulated).shape)\n",
    "print('n_subjects, n_TRs, n_voxels')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "1cf43252",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The data shape you given, this is check the data shape whether it is right or not\n",
      "(803, 360180, 4)\n"
     ]
    }
   ],
   "source": [
    "print('The data shape you given, this is check the data shape whether it is right or not')\n",
    "print(data.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "5e6bcb60",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(4, 803, 360180)\n"
     ]
    }
   ],
   "source": [
    "data_trans = np.transpose(data, (2,0,1))\n",
    "print(data_trans.shape)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d9a6a67c",
   "metadata": {},
   "source": [
    "## ISC movie comparison"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9b6d5750",
   "metadata": {},
   "outputs": [],
   "source": [
    "iscs = isc(data_trans, pairwise=True, summary_statistic=None)\n",
    "observed, p, distribution = permutation_isc(iscs, pairwise=True,\n",
    "                                            summary_statistic='median',\n",
    "                                            n_permutations=200)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5f0ed4ba",
   "metadata": {},
   "outputs": [],
   "source": [
    "print(f\"Null distribution shape = {distribution.shape}\"\n",
    "      f\"\\ni.e., {distribution.shape[0]} permutations \"\n",
    "      f\"and {distribution.shape[1]} voxels\")\n",
    "\n",
    "# Get actual ISC value and p-value for first voxel\n",
    "print(f\"Actual observed ISC value for first voxel = {observed[0][0]:.3f},\"\n",
    "      f\"\\np-value from permutation test = {p[0]:.3f}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3ce8d816",
   "metadata": {},
   "source": [
    "## group comparison\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "51df2819",
   "metadata": {},
   "outputs": [],
   "source": [
    "group1 = [data_trans[0],data_trans[1]]\n",
    "group2 = [data_trans[2],data_trans[3]]\n",
    "group_assignment = group1 + group2 \n",
    "iscs = isc(data, pairwise=True, summary_statistic=None)\n",
    "observed, p, distribution = permutation_isc(iscs,\n",
    "                                            group_assignment=group_assignment,\n",
    "                                            pairwise=True,\n",
    "                                            summary_statistic='median',\n",
    "                                            n_permutations=4697184)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "12dab581",
   "metadata": {},
   "outputs": [],
   "source": [
    "print(f\"Null distribution shape = {distribution.shape}\"\n",
    "      f\"\\ni.e., {distribution.shape[0]} permutations \"\n",
    "      f\"and {distribution.shape[1]} voxels\")\n",
    "\n",
    "# Get actual ISC value and p-value for first voxel\n",
    "print(f\"Actual observed group difference in ISC values \"\n",
    "      f\"for first voxel = {observed[0]:.3f},\"\n",
    "      f\"\\np-value from permutation test = {p[0]:.3f}\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c2c30fbf",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c295ef1c",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "284dd15b",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python [conda env:Neuro] *",
   "language": "python",
   "name": "conda-env-Neuro-py"
  },
  "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.13"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}