2_analyze_phenotypic_data.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "26b163e1",
   "metadata": {},
   "source": [
    "# Analyze phenotypic data\n",
    "### Start by importing packages and loading the dataframe from the previous step"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "eaef0543",
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import pingouin\n",
    "import seaborn as sns\n",
    "\n",
    "df = pd.read_pickle('df.pkl')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "29be3f90",
   "metadata": {},
   "source": [
    "### Check for group differences (Table 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6a533b38",
   "metadata": {},
   "outputs": [],
   "source": [
    "# PRINT COHORT-WIDE STATISTICS\n",
    "print('AVERAGE AGE:',np.mean(df['AGE']),'(',np.std(df['AGE']),')')\n",
    "print('AVERAGE HAND:',np.mean(df['HAND']),'(',np.std(df['HAND']),')')\n",
    "print('AVERAGE gFA:',np.mean(df['gFA']),'(',np.std(df['gFA']),')')\n",
    "print('AVERAGE TBV:',np.mean(df['ICV']),'(',np.std(df['ICV']),') \\n')\n",
    "\n",
    "# Get the mean, sum, and std of stats broken down by groups\n",
    "mean_stats_ht_lt = df.groupby('TR_RD').mean()\n",
    "sum_stats_ht_lt = df.groupby('TR_RD').sum()\n",
    "sd_stats_ht_lt = df.groupby('TR_RD').std()\n",
    "\n",
    "# Get indices of high and low TOWRE participants\n",
    "inds_low_t = (df['TR_RD']==1)\n",
    "inds_high_t = (df['TR_RD']==0)\n",
    "\n",
    "# Get number of male and females in the low and high TOWRE groups\n",
    "inds_rd = (df['TR_RD']==1)\n",
    "inds_tr = (df['TR_RD']==0)\n",
    "m_rd = sum(inds_rd) - sum(df[inds_rd]['SEX'])\n",
    "f_rd = sum(inds_rd) - m_rd\n",
    "m_tr = sum(inds_tr) - sum(df[inds_tr]['SEX'])\n",
    "f_tr = sum(inds_tr) - m_tr\n",
    "print(m_rd, f_rd, m_tr, f_tr)\n",
    "\n",
    "# DIFFERENCES IN SEX (Chi-Square)\n",
    "#print('SEX: M HIGH T:',m_high_t,', F HIGH T:',f_high_t,', M LOW T',m_low_t,', F LOW T',f_low_t,'\\n')\n",
    "\n",
    "# DIFFERENCES IN AGE (Welch's t-test)\n",
    "print('COMPARE AGE: HIGH T MEAN (SD):',mean_stats_ht_lt['AGE'][0],'(',sd_stats_ht_lt['AGE'][0],'), LOW T MEAN (SD):', \n",
    "      mean_stats_ht_lt['AGE'][1],'(',sd_stats_ht_lt['AGE'][1],')\\n')\n",
    "print(str(pingouin.ttest(df['AGE'][inds_high_t], df['AGE'][inds_low_t]))+'\\n')\n",
    "\n",
    "# DIFFERENCES IN HANDEDNESS (Welch's t-test)\n",
    "print('COMPARE HAND: HIGH T MEAN (SD):',mean_stats_ht_lt['HAND'][0],'(',sd_stats_ht_lt['HAND'][0],'), LOW T MEAN (SD):', \n",
    "      mean_stats_ht_lt['HAND'][1],'(',sd_stats_ht_lt['HAND'][1],')\\n')\n",
    "print(str(pingouin.ttest(df['HAND'][inds_high_t], df['HAND'][inds_low_t]))+'\\n')\n",
    "\n",
    "# DIFFERENCES IN gFA (Welch's t-test)\n",
    "print('COMPARE gFA: HIGH T MEAN (SD):',mean_stats_ht_lt['gFA'][0],'(',sd_stats_ht_lt['gFA'][0],'), LOW T MEAN (SD):', \n",
    "      mean_stats_ht_lt['gFA'][1],'(',sd_stats_ht_lt['gFA'][1],')\\n')\n",
    "print(str(pingouin.ttest(df['gFA'][inds_high_t], df['gFA'][inds_low_t]))+'\\n')\n",
    "\n",
    "# DIFFERENCES IN ICV (Welch's t-test)\n",
    "print('COMPARE ICV: HIGH T MEAN (SD):',mean_stats_ht_lt['ICV'][0],'(',sd_stats_ht_lt['ICV'][0],'), LOW T MEAN (SD):', \n",
    "      mean_stats_ht_lt['ICV'][1],'(',sd_stats_ht_lt['ICV'][1],')\\n')\n",
    "print(str(pingouin.ttest(df['ICV'][inds_high_t], df['ICV'][inds_low_t])))\n",
    "\n",
    "# DIFFERENCES IN TBV (Welch's t-test)\n",
    "print('COMPARE TBV: HIGH T MEAN (SD):',mean_stats_ht_lt['TBV'][0],'(',sd_stats_ht_lt['TBV'][0],'), LOW T MEAN (SD):', \n",
    "      mean_stats_ht_lt['TBV'][1],'(',sd_stats_ht_lt['TBV'][1],')\\n')\n",
    "print(str(pingouin.ttest(df['TBV'][inds_high_t], df['TBV'][inds_low_t])))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2098c068",
   "metadata": {},
   "source": [
    "### Plot TOWRE Distribution (Figure 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "20b2fb7b",
   "metadata": {},
   "outputs": [],
   "source": [
    "font = {'family' : 'Arial',\n",
    "        'weight' : 'bold',\n",
    "        'size'   : 22}\n",
    "sns.set(font_scale=1.3)\n",
    "sns.set_style(\"white\")\n",
    "plt.rc('font', **font)\n",
    "t_plot = sns.displot(df, x=\"TOWRE\", hue=\"RD_TR\",bins=20, palette=['teal','red'],multiple=\"stack\")\n",
    "t_plot._legend.remove()\n",
    "plt.vlines(85,0,90,'k',linewidth=5,linestyle=':')\n",
    "plt.ylim([0,85])\n",
    "plt.xlabel('TOWRE Total Score, Standardized',fontsize=20,fontweight='bold')\n",
    "plt.ylabel('Number of Participants',fontsize=20,fontweight='bold')\n",
    "plt.savefig('TOWRE.pdf') # use Adobe Illustrator to convert to EPS in order to keep transparent background\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cf3ec6be-d305-49a1-9cf0-094c257ec475",
   "metadata": {},
   "source": [
    "## Here is where you can see how different properties covary with each other (Figure 2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "58983ec9-1957-4a8e-bb73-0ecac18f9923",
   "metadata": {},
   "outputs": [],
   "source": [
    "sns.set(font_scale = 2)\n",
    "sns.set_style(\"white\") # Format background style\n",
    "# Get dataframe with just covariates\n",
    "df_covars = df[['AGE','SEX','HAND','gFA','ICV','TBV','TOWRE']]\n",
    "\n",
    "# Calculate correlations across all pairwise columns\n",
    "corr = df_covars.corr(method='spearman')\n",
    "\n",
    "# Generate a mask for the upper triangle\n",
    "mask = np.triu(np.ones_like(corr, dtype=bool))\n",
    "\n",
    "# Set up the matplotlib figure\n",
    "f, ax = plt.subplots(figsize=(11, 9))\n",
    "\n",
    "# Generate a custom diverging colormap\n",
    "cmap = sns.diverging_palette(230, 20, as_cmap=True)\n",
    "\n",
    "# Draw the heatmap with the mask and correct aspect ratio\n",
    "sns.heatmap(corr, mask=mask, cmap=cmap, vmax = .45, vmin = -.45, center=0,\n",
    "            square=True, linewidths=.5, cbar_kws={\"shrink\": .5}, annot=True, annot_kws={\"fontsize\":15})\n",
    "#plt.savefig('covars_corr.eps',format='eps')\n",
    "plt.show()\n",
    "pingouin.rcorr(df_covars,method='spearman',padjust='fdr_bh')"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "nipype",
   "language": "python",
   "name": "nipype"
  },
  "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.8.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}