FINAL_DGN_OLS.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "import statsmodels.api as sm\n",
    "from statsmodels.stats.multitest import fdrcorrection"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>DQA1*01:01</th>\n",
       "      <th>DQA1*01:02</th>\n",
       "      <th>DQA1*01:03</th>\n",
       "      <th>DQA1*01:04</th>\n",
       "      <th>DQA1*01:05</th>\n",
       "      <th>DQA1*02:01</th>\n",
       "      <th>DQA1*03:01</th>\n",
       "      <th>DQA1*03:02</th>\n",
       "      <th>DQA1*03:03</th>\n",
       "      <th>DQA1*04:01</th>\n",
       "      <th>...</th>\n",
       "      <th>DQB1*04:02</th>\n",
       "      <th>DQB1*05:01</th>\n",
       "      <th>DQB1*05:02</th>\n",
       "      <th>DQB1*05:03</th>\n",
       "      <th>DQB1*05:04</th>\n",
       "      <th>DQB1*06:01</th>\n",
       "      <th>DQB1*06:02</th>\n",
       "      <th>DQB1*06:03</th>\n",
       "      <th>DQB1*06:04</th>\n",
       "      <th>DQB1*06:09</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>patid</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>LD0014</th>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD0041</th>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD0038</th>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD0084</th>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD0022</th>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>...</th>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD1282</th>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD1271</th>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD1252</th>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD0165</th>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD0102</th>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>895 rows × 32 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "        DQA1*01:01  DQA1*01:02  DQA1*01:03  DQA1*01:04  DQA1*01:05  \\\n",
       "patid                                                                \n",
       "LD0014         1.0         0.0         0.0         0.0         0.0   \n",
       "LD0041         0.0         1.0         0.0         0.0         0.0   \n",
       "LD0038         0.0         1.0         0.0         0.0         0.0   \n",
       "LD0084         0.0         1.0         1.0         0.0         0.0   \n",
       "LD0022         1.0         0.0         0.0         0.0         0.0   \n",
       "...            ...         ...         ...         ...         ...   \n",
       "LD1282         0.0         0.0         0.0         0.0         0.0   \n",
       "LD1271         0.0         0.0         1.0         0.0         0.0   \n",
       "LD1252         0.0         1.0         1.0         0.0         0.0   \n",
       "LD0165         0.0         1.0         0.0         0.0         0.0   \n",
       "LD0102         0.0         1.0         0.0         0.0         0.0   \n",
       "\n",
       "        DQA1*02:01  DQA1*03:01  DQA1*03:02  DQA1*03:03  DQA1*04:01  ...  \\\n",
       "patid                                                               ...   \n",
       "LD0014         0.0         1.0         0.0         0.0         0.0  ...   \n",
       "LD0041         0.0         1.0         0.0         0.0         0.0  ...   \n",
       "LD0038         0.0         0.0         0.0         0.0         0.0  ...   \n",
       "LD0084         0.0         0.0         0.0         0.0         0.0  ...   \n",
       "LD0022         0.0         1.0         0.0         0.0         0.0  ...   \n",
       "...            ...         ...         ...         ...         ...  ...   \n",
       "LD1282         0.0         0.0         0.0         1.0         0.0  ...   \n",
       "LD1271         0.0         0.0         1.0         0.0         0.0  ...   \n",
       "LD1252         0.0         0.0         0.0         0.0         0.0  ...   \n",
       "LD0165         0.0         1.0         0.0         0.0         0.0  ...   \n",
       "LD0102         0.0         0.0         0.0         0.0         0.0  ...   \n",
       "\n",
       "        DQB1*04:02  DQB1*05:01  DQB1*05:02  DQB1*05:03  DQB1*05:04  \\\n",
       "patid                                                                \n",
       "LD0014         0.0         1.0         0.0         0.0         0.0   \n",
       "LD0041         0.0         0.0         0.0         0.0         0.0   \n",
       "LD0038         0.0         0.0         0.0         0.0         0.0   \n",
       "LD0084         0.0         0.0         0.0         0.0         0.0   \n",
       "LD0022         0.0         1.0         0.0         0.0         0.0   \n",
       "...            ...         ...         ...         ...         ...   \n",
       "LD1282         0.0         0.0         0.0         0.0         0.0   \n",
       "LD1271         0.0         0.0         0.0         0.0         0.0   \n",
       "LD1252         0.0         0.0         0.0         0.0         0.0   \n",
       "LD0165         0.0         0.0         0.0         0.0         0.0   \n",
       "LD0102         0.0         0.0         0.0         0.0         0.0   \n",
       "\n",
       "        DQB1*06:01  DQB1*06:02  DQB1*06:03  DQB1*06:04  DQB1*06:09  \n",
       "patid                                                               \n",
       "LD0014         0.0         0.0         0.0         0.0         0.0  \n",
       "LD0041         0.0         1.0         0.0         0.0         0.0  \n",
       "LD0038         0.0         1.0         0.0         0.0         0.0  \n",
       "LD0084         1.0         1.0         0.0         0.0         0.0  \n",
       "LD0022         0.0         0.0         0.0         0.0         0.0  \n",
       "...            ...         ...         ...         ...         ...  \n",
       "LD1282         0.0         0.0         0.0         0.0         0.0  \n",
       "LD1271         0.0         0.0         1.0         0.0         0.0  \n",
       "LD1252         0.0         1.0         1.0         0.0         0.0  \n",
       "LD0165         0.0         0.0         0.0         1.0         0.0  \n",
       "LD0102         0.0         1.0         0.0         0.0         0.0  \n",
       "\n",
       "[895 rows x 32 columns]"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "HLA_table = pd.read_csv(\"DGN_HLA_df.csv\").set_index(\"patid\").apply(lambda x: x > 0).astype(\"float\")\n",
    "HLA_table"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>TRAV1-1</th>\n",
       "      <th>TRAV1-2</th>\n",
       "      <th>TRAV10</th>\n",
       "      <th>TRAV12-1</th>\n",
       "      <th>TRAV12-2</th>\n",
       "      <th>TRAV12-3</th>\n",
       "      <th>TRAV13-1</th>\n",
       "      <th>TRAV13-2</th>\n",
       "      <th>TRAV14/DV4</th>\n",
       "      <th>TRAV16</th>\n",
       "      <th>...</th>\n",
       "      <th>TRAV5</th>\n",
       "      <th>TRAV6</th>\n",
       "      <th>TRAV7</th>\n",
       "      <th>TRAV8-1</th>\n",
       "      <th>TRAV8-2</th>\n",
       "      <th>TRAV8-3</th>\n",
       "      <th>TRAV8-4</th>\n",
       "      <th>TRAV8-6</th>\n",
       "      <th>TRAV9-1</th>\n",
       "      <th>TRAV9-2</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>patid</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>LD0001</th>\n",
       "      <td>60.0</td>\n",
       "      <td>125.0</td>\n",
       "      <td>94.0</td>\n",
       "      <td>220.0</td>\n",
       "      <td>371.0</td>\n",
       "      <td>202.0</td>\n",
       "      <td>333.0</td>\n",
       "      <td>132.0</td>\n",
       "      <td>187.0</td>\n",
       "      <td>90.0</td>\n",
       "      <td>...</td>\n",
       "      <td>79.0</td>\n",
       "      <td>78.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>107.0</td>\n",
       "      <td>60.0</td>\n",
       "      <td>169.0</td>\n",
       "      <td>241.0</td>\n",
       "      <td>157.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>409.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD0002</th>\n",
       "      <td>66.0</td>\n",
       "      <td>217.0</td>\n",
       "      <td>113.0</td>\n",
       "      <td>409.0</td>\n",
       "      <td>376.0</td>\n",
       "      <td>243.0</td>\n",
       "      <td>500.0</td>\n",
       "      <td>201.0</td>\n",
       "      <td>187.0</td>\n",
       "      <td>119.0</td>\n",
       "      <td>...</td>\n",
       "      <td>76.0</td>\n",
       "      <td>76.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>96.0</td>\n",
       "      <td>90.0</td>\n",
       "      <td>213.0</td>\n",
       "      <td>313.0</td>\n",
       "      <td>240.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>527.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD0003</th>\n",
       "      <td>118.0</td>\n",
       "      <td>207.0</td>\n",
       "      <td>82.0</td>\n",
       "      <td>238.0</td>\n",
       "      <td>571.0</td>\n",
       "      <td>269.0</td>\n",
       "      <td>437.0</td>\n",
       "      <td>135.0</td>\n",
       "      <td>263.0</td>\n",
       "      <td>103.0</td>\n",
       "      <td>...</td>\n",
       "      <td>103.0</td>\n",
       "      <td>61.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>73.0</td>\n",
       "      <td>101.0</td>\n",
       "      <td>153.0</td>\n",
       "      <td>280.0</td>\n",
       "      <td>185.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>330.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD0006</th>\n",
       "      <td>106.0</td>\n",
       "      <td>256.0</td>\n",
       "      <td>110.0</td>\n",
       "      <td>304.0</td>\n",
       "      <td>312.0</td>\n",
       "      <td>246.0</td>\n",
       "      <td>439.0</td>\n",
       "      <td>209.0</td>\n",
       "      <td>169.0</td>\n",
       "      <td>147.0</td>\n",
       "      <td>...</td>\n",
       "      <td>96.0</td>\n",
       "      <td>133.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>120.0</td>\n",
       "      <td>99.0</td>\n",
       "      <td>194.0</td>\n",
       "      <td>294.0</td>\n",
       "      <td>214.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>412.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD0007</th>\n",
       "      <td>75.0</td>\n",
       "      <td>202.0</td>\n",
       "      <td>114.0</td>\n",
       "      <td>408.0</td>\n",
       "      <td>352.0</td>\n",
       "      <td>276.0</td>\n",
       "      <td>693.0</td>\n",
       "      <td>280.0</td>\n",
       "      <td>364.0</td>\n",
       "      <td>138.0</td>\n",
       "      <td>...</td>\n",
       "      <td>120.0</td>\n",
       "      <td>129.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>114.0</td>\n",
       "      <td>110.0</td>\n",
       "      <td>244.0</td>\n",
       "      <td>360.0</td>\n",
       "      <td>182.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>689.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>...</th>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD1357</th>\n",
       "      <td>89.0</td>\n",
       "      <td>260.0</td>\n",
       "      <td>80.0</td>\n",
       "      <td>191.0</td>\n",
       "      <td>287.0</td>\n",
       "      <td>164.0</td>\n",
       "      <td>281.0</td>\n",
       "      <td>121.0</td>\n",
       "      <td>98.0</td>\n",
       "      <td>69.0</td>\n",
       "      <td>...</td>\n",
       "      <td>70.0</td>\n",
       "      <td>39.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>71.0</td>\n",
       "      <td>57.0</td>\n",
       "      <td>101.0</td>\n",
       "      <td>194.0</td>\n",
       "      <td>167.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>295.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD1361</th>\n",
       "      <td>66.0</td>\n",
       "      <td>131.0</td>\n",
       "      <td>170.0</td>\n",
       "      <td>329.0</td>\n",
       "      <td>291.0</td>\n",
       "      <td>257.0</td>\n",
       "      <td>650.0</td>\n",
       "      <td>346.0</td>\n",
       "      <td>257.0</td>\n",
       "      <td>164.0</td>\n",
       "      <td>...</td>\n",
       "      <td>123.0</td>\n",
       "      <td>114.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>79.0</td>\n",
       "      <td>108.0</td>\n",
       "      <td>217.0</td>\n",
       "      <td>307.0</td>\n",
       "      <td>268.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>581.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD1362</th>\n",
       "      <td>59.0</td>\n",
       "      <td>293.0</td>\n",
       "      <td>52.0</td>\n",
       "      <td>265.0</td>\n",
       "      <td>256.0</td>\n",
       "      <td>205.0</td>\n",
       "      <td>361.0</td>\n",
       "      <td>136.0</td>\n",
       "      <td>155.0</td>\n",
       "      <td>102.0</td>\n",
       "      <td>...</td>\n",
       "      <td>111.0</td>\n",
       "      <td>79.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>63.0</td>\n",
       "      <td>57.0</td>\n",
       "      <td>142.0</td>\n",
       "      <td>242.0</td>\n",
       "      <td>159.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>238.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD1364</th>\n",
       "      <td>72.0</td>\n",
       "      <td>152.0</td>\n",
       "      <td>123.0</td>\n",
       "      <td>196.0</td>\n",
       "      <td>346.0</td>\n",
       "      <td>200.0</td>\n",
       "      <td>331.0</td>\n",
       "      <td>138.0</td>\n",
       "      <td>164.0</td>\n",
       "      <td>97.0</td>\n",
       "      <td>...</td>\n",
       "      <td>100.0</td>\n",
       "      <td>97.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>72.0</td>\n",
       "      <td>48.0</td>\n",
       "      <td>151.0</td>\n",
       "      <td>225.0</td>\n",
       "      <td>158.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>280.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD1366</th>\n",
       "      <td>135.0</td>\n",
       "      <td>281.0</td>\n",
       "      <td>80.0</td>\n",
       "      <td>289.0</td>\n",
       "      <td>251.0</td>\n",
       "      <td>343.0</td>\n",
       "      <td>431.0</td>\n",
       "      <td>179.0</td>\n",
       "      <td>280.0</td>\n",
       "      <td>132.0</td>\n",
       "      <td>...</td>\n",
       "      <td>104.0</td>\n",
       "      <td>87.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>110.0</td>\n",
       "      <td>96.0</td>\n",
       "      <td>342.0</td>\n",
       "      <td>246.0</td>\n",
       "      <td>113.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>295.0</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>895 rows × 45 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "        TRAV1-1  TRAV1-2  TRAV10  TRAV12-1  TRAV12-2  TRAV12-3  TRAV13-1  \\\n",
       "patid                                                                      \n",
       "LD0001     60.0    125.0    94.0     220.0     371.0     202.0     333.0   \n",
       "LD0002     66.0    217.0   113.0     409.0     376.0     243.0     500.0   \n",
       "LD0003    118.0    207.0    82.0     238.0     571.0     269.0     437.0   \n",
       "LD0006    106.0    256.0   110.0     304.0     312.0     246.0     439.0   \n",
       "LD0007     75.0    202.0   114.0     408.0     352.0     276.0     693.0   \n",
       "...         ...      ...     ...       ...       ...       ...       ...   \n",
       "LD1357     89.0    260.0    80.0     191.0     287.0     164.0     281.0   \n",
       "LD1361     66.0    131.0   170.0     329.0     291.0     257.0     650.0   \n",
       "LD1362     59.0    293.0    52.0     265.0     256.0     205.0     361.0   \n",
       "LD1364     72.0    152.0   123.0     196.0     346.0     200.0     331.0   \n",
       "LD1366    135.0    281.0    80.0     289.0     251.0     343.0     431.0   \n",
       "\n",
       "        TRAV13-2  TRAV14/DV4  TRAV16  ...  TRAV5  TRAV6  TRAV7  TRAV8-1  \\\n",
       "patid                                 ...                                 \n",
       "LD0001     132.0       187.0    90.0  ...   79.0   78.0    0.0    107.0   \n",
       "LD0002     201.0       187.0   119.0  ...   76.0   76.0    0.0     96.0   \n",
       "LD0003     135.0       263.0   103.0  ...  103.0   61.0    0.0     73.0   \n",
       "LD0006     209.0       169.0   147.0  ...   96.0  133.0    0.0    120.0   \n",
       "LD0007     280.0       364.0   138.0  ...  120.0  129.0    0.0    114.0   \n",
       "...          ...         ...     ...  ...    ...    ...    ...      ...   \n",
       "LD1357     121.0        98.0    69.0  ...   70.0   39.0    0.0     71.0   \n",
       "LD1361     346.0       257.0   164.0  ...  123.0  114.0    0.0     79.0   \n",
       "LD1362     136.0       155.0   102.0  ...  111.0   79.0    0.0     63.0   \n",
       "LD1364     138.0       164.0    97.0  ...  100.0   97.0    1.0     72.0   \n",
       "LD1366     179.0       280.0   132.0  ...  104.0   87.0    1.0    110.0   \n",
       "\n",
       "        TRAV8-2  TRAV8-3  TRAV8-4  TRAV8-6  TRAV9-1  TRAV9-2  \n",
       "patid                                                         \n",
       "LD0001     60.0    169.0    241.0    157.0      0.0    409.0  \n",
       "LD0002     90.0    213.0    313.0    240.0      1.0    527.0  \n",
       "LD0003    101.0    153.0    280.0    185.0      0.0    330.0  \n",
       "LD0006     99.0    194.0    294.0    214.0      1.0    412.0  \n",
       "LD0007    110.0    244.0    360.0    182.0      0.0    689.0  \n",
       "...         ...      ...      ...      ...      ...      ...  \n",
       "LD1357     57.0    101.0    194.0    167.0      0.0    295.0  \n",
       "LD1361    108.0    217.0    307.0    268.0      0.0    581.0  \n",
       "LD1362     57.0    142.0    242.0    159.0      0.0    238.0  \n",
       "LD1364     48.0    151.0    225.0    158.0      0.0    280.0  \n",
       "LD1366     96.0    342.0    246.0    113.0      0.0    295.0  \n",
       "\n",
       "[895 rows x 45 columns]"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "family_counts_table = pd.read_csv(\"DGN_family_counts_table.csv\")\n",
    "family_counts_table = family_counts_table.set_index(\"patid\")\n",
    "family_counts_table = family_counts_table[[x for x in family_counts_table.columns if x.startswith(\"TRAV\")]]\n",
    "family_counts_table"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "tcr_columns = list(family_counts_table.columns)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "CDR_NUM_AA = 22 # 22 positions, first 12 are CDR1, then 10 are CDR2\n",
    "DGN_NUM_PATIENTS = 895"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "TRV = \"TRBV\" # TRAV\n",
    "\n",
    "positions_AA_count_dfs = []\n",
    "for i in range(CDR_NUM_AA):\n",
    "    positions_AA_count_dfs.append(pd.read_csv(\"amino_acids_correlation/DGN_{}_CDR1_CDR2_position_{}_AA_counts_table.csv\".format(TRV, i)).set_index(\"patid\"))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>D</th>\n",
       "      <th>E</th>\n",
       "      <th>F</th>\n",
       "      <th>G</th>\n",
       "      <th>N</th>\n",
       "      <th>Q</th>\n",
       "      <th>R</th>\n",
       "      <th>S</th>\n",
       "      <th>T</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>patid</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>LD0001</th>\n",
       "      <td>253.0</td>\n",
       "      <td>502.0</td>\n",
       "      <td>915.0</td>\n",
       "      <td>4322.0</td>\n",
       "      <td>2063.0</td>\n",
       "      <td>498.0</td>\n",
       "      <td>105.0</td>\n",
       "      <td>81.0</td>\n",
       "      <td>89.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD0002</th>\n",
       "      <td>480.0</td>\n",
       "      <td>807.0</td>\n",
       "      <td>1035.0</td>\n",
       "      <td>5966.0</td>\n",
       "      <td>3591.0</td>\n",
       "      <td>747.0</td>\n",
       "      <td>224.0</td>\n",
       "      <td>163.0</td>\n",
       "      <td>173.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD0003</th>\n",
       "      <td>144.0</td>\n",
       "      <td>644.0</td>\n",
       "      <td>594.0</td>\n",
       "      <td>4553.0</td>\n",
       "      <td>2780.0</td>\n",
       "      <td>307.0</td>\n",
       "      <td>112.0</td>\n",
       "      <td>119.0</td>\n",
       "      <td>256.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD0006</th>\n",
       "      <td>739.0</td>\n",
       "      <td>687.0</td>\n",
       "      <td>1416.0</td>\n",
       "      <td>5127.0</td>\n",
       "      <td>3195.0</td>\n",
       "      <td>815.0</td>\n",
       "      <td>194.0</td>\n",
       "      <td>107.0</td>\n",
       "      <td>163.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD0007</th>\n",
       "      <td>606.0</td>\n",
       "      <td>891.0</td>\n",
       "      <td>1298.0</td>\n",
       "      <td>7370.0</td>\n",
       "      <td>3882.0</td>\n",
       "      <td>758.0</td>\n",
       "      <td>138.0</td>\n",
       "      <td>138.0</td>\n",
       "      <td>220.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>...</th>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD1357</th>\n",
       "      <td>331.0</td>\n",
       "      <td>445.0</td>\n",
       "      <td>760.0</td>\n",
       "      <td>4458.0</td>\n",
       "      <td>2111.0</td>\n",
       "      <td>654.0</td>\n",
       "      <td>185.0</td>\n",
       "      <td>55.0</td>\n",
       "      <td>287.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD1361</th>\n",
       "      <td>277.0</td>\n",
       "      <td>1075.0</td>\n",
       "      <td>1408.0</td>\n",
       "      <td>7654.0</td>\n",
       "      <td>4325.0</td>\n",
       "      <td>1189.0</td>\n",
       "      <td>183.0</td>\n",
       "      <td>186.0</td>\n",
       "      <td>227.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD1362</th>\n",
       "      <td>613.0</td>\n",
       "      <td>731.0</td>\n",
       "      <td>764.0</td>\n",
       "      <td>4410.0</td>\n",
       "      <td>2579.0</td>\n",
       "      <td>511.0</td>\n",
       "      <td>126.0</td>\n",
       "      <td>76.0</td>\n",
       "      <td>216.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD1364</th>\n",
       "      <td>90.0</td>\n",
       "      <td>839.0</td>\n",
       "      <td>743.0</td>\n",
       "      <td>6670.0</td>\n",
       "      <td>3625.0</td>\n",
       "      <td>570.0</td>\n",
       "      <td>124.0</td>\n",
       "      <td>209.0</td>\n",
       "      <td>387.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD1366</th>\n",
       "      <td>681.0</td>\n",
       "      <td>1079.0</td>\n",
       "      <td>1292.0</td>\n",
       "      <td>8514.0</td>\n",
       "      <td>3533.0</td>\n",
       "      <td>758.0</td>\n",
       "      <td>186.0</td>\n",
       "      <td>74.0</td>\n",
       "      <td>445.0</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>895 rows × 9 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "            D       E       F       G       N       Q      R      S      T\n",
       "patid                                                                     \n",
       "LD0001  253.0   502.0   915.0  4322.0  2063.0   498.0  105.0   81.0   89.0\n",
       "LD0002  480.0   807.0  1035.0  5966.0  3591.0   747.0  224.0  163.0  173.0\n",
       "LD0003  144.0   644.0   594.0  4553.0  2780.0   307.0  112.0  119.0  256.0\n",
       "LD0006  739.0   687.0  1416.0  5127.0  3195.0   815.0  194.0  107.0  163.0\n",
       "LD0007  606.0   891.0  1298.0  7370.0  3882.0   758.0  138.0  138.0  220.0\n",
       "...       ...     ...     ...     ...     ...     ...    ...    ...    ...\n",
       "LD1357  331.0   445.0   760.0  4458.0  2111.0   654.0  185.0   55.0  287.0\n",
       "LD1361  277.0  1075.0  1408.0  7654.0  4325.0  1189.0  183.0  186.0  227.0\n",
       "LD1362  613.0   731.0   764.0  4410.0  2579.0   511.0  126.0   76.0  216.0\n",
       "LD1364   90.0   839.0   743.0  6670.0  3625.0   570.0  124.0  209.0  387.0\n",
       "LD1366  681.0  1079.0  1292.0  8514.0  3533.0   758.0  186.0   74.0  445.0\n",
       "\n",
       "[895 rows x 9 columns]"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "positions_AA_count_dfs[1]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Using family counts and TRAV sequences, for every patient, get count of each AA at each position of CDR2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "def counts_to_usage_table(counts_table, drop_columns={}):\n",
    "    counts_table = counts_table.drop(axis=1, labels=drop_columns)\n",
    "    usage_table = counts_table.div(counts_table.sum(axis=1), axis=0)\n",
    "    return usage_table"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "family_usage_table = counts_to_usage_table(family_counts_table)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>TRAV1-1</th>\n",
       "      <th>TRAV1-2</th>\n",
       "      <th>TRAV10</th>\n",
       "      <th>TRAV12-1</th>\n",
       "      <th>TRAV12-2</th>\n",
       "      <th>TRAV12-3</th>\n",
       "      <th>TRAV13-1</th>\n",
       "      <th>TRAV13-2</th>\n",
       "      <th>TRAV14/DV4</th>\n",
       "      <th>TRAV16</th>\n",
       "      <th>...</th>\n",
       "      <th>TRAV5</th>\n",
       "      <th>TRAV6</th>\n",
       "      <th>TRAV7</th>\n",
       "      <th>TRAV8-1</th>\n",
       "      <th>TRAV8-2</th>\n",
       "      <th>TRAV8-3</th>\n",
       "      <th>TRAV8-4</th>\n",
       "      <th>TRAV8-6</th>\n",
       "      <th>TRAV9-1</th>\n",
       "      <th>TRAV9-2</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>patid</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>LD0001</th>\n",
       "      <td>0.010435</td>\n",
       "      <td>0.021739</td>\n",
       "      <td>0.016348</td>\n",
       "      <td>0.038261</td>\n",
       "      <td>0.064522</td>\n",
       "      <td>0.035130</td>\n",
       "      <td>0.057913</td>\n",
       "      <td>0.022957</td>\n",
       "      <td>0.032522</td>\n",
       "      <td>0.015652</td>\n",
       "      <td>...</td>\n",
       "      <td>0.013739</td>\n",
       "      <td>0.013565</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.018609</td>\n",
       "      <td>0.010435</td>\n",
       "      <td>0.029391</td>\n",
       "      <td>0.041913</td>\n",
       "      <td>0.027304</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.071130</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD0002</th>\n",
       "      <td>0.007591</td>\n",
       "      <td>0.024957</td>\n",
       "      <td>0.012996</td>\n",
       "      <td>0.047039</td>\n",
       "      <td>0.043243</td>\n",
       "      <td>0.027947</td>\n",
       "      <td>0.057504</td>\n",
       "      <td>0.023117</td>\n",
       "      <td>0.021507</td>\n",
       "      <td>0.013686</td>\n",
       "      <td>...</td>\n",
       "      <td>0.008741</td>\n",
       "      <td>0.008741</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.011041</td>\n",
       "      <td>0.010351</td>\n",
       "      <td>0.024497</td>\n",
       "      <td>0.035998</td>\n",
       "      <td>0.027602</td>\n",
       "      <td>0.000115</td>\n",
       "      <td>0.060610</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD0003</th>\n",
       "      <td>0.016874</td>\n",
       "      <td>0.029601</td>\n",
       "      <td>0.011726</td>\n",
       "      <td>0.034034</td>\n",
       "      <td>0.081653</td>\n",
       "      <td>0.038467</td>\n",
       "      <td>0.062491</td>\n",
       "      <td>0.019305</td>\n",
       "      <td>0.037609</td>\n",
       "      <td>0.014729</td>\n",
       "      <td>...</td>\n",
       "      <td>0.014729</td>\n",
       "      <td>0.008723</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.010439</td>\n",
       "      <td>0.014443</td>\n",
       "      <td>0.021879</td>\n",
       "      <td>0.040040</td>\n",
       "      <td>0.026455</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.047190</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD0006</th>\n",
       "      <td>0.013472</td>\n",
       "      <td>0.032537</td>\n",
       "      <td>0.013981</td>\n",
       "      <td>0.038638</td>\n",
       "      <td>0.039654</td>\n",
       "      <td>0.031266</td>\n",
       "      <td>0.055796</td>\n",
       "      <td>0.026563</td>\n",
       "      <td>0.021479</td>\n",
       "      <td>0.018683</td>\n",
       "      <td>...</td>\n",
       "      <td>0.012201</td>\n",
       "      <td>0.016904</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.015252</td>\n",
       "      <td>0.012583</td>\n",
       "      <td>0.024657</td>\n",
       "      <td>0.037367</td>\n",
       "      <td>0.027199</td>\n",
       "      <td>0.000127</td>\n",
       "      <td>0.052364</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD0007</th>\n",
       "      <td>0.008087</td>\n",
       "      <td>0.021781</td>\n",
       "      <td>0.012292</td>\n",
       "      <td>0.043994</td>\n",
       "      <td>0.037956</td>\n",
       "      <td>0.029761</td>\n",
       "      <td>0.074725</td>\n",
       "      <td>0.030192</td>\n",
       "      <td>0.039250</td>\n",
       "      <td>0.014880</td>\n",
       "      <td>...</td>\n",
       "      <td>0.012939</td>\n",
       "      <td>0.013910</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.012292</td>\n",
       "      <td>0.011861</td>\n",
       "      <td>0.026310</td>\n",
       "      <td>0.038818</td>\n",
       "      <td>0.019625</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.074294</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>...</th>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD1357</th>\n",
       "      <td>0.017994</td>\n",
       "      <td>0.052568</td>\n",
       "      <td>0.016175</td>\n",
       "      <td>0.038617</td>\n",
       "      <td>0.058027</td>\n",
       "      <td>0.033158</td>\n",
       "      <td>0.056814</td>\n",
       "      <td>0.024464</td>\n",
       "      <td>0.019814</td>\n",
       "      <td>0.013951</td>\n",
       "      <td>...</td>\n",
       "      <td>0.014153</td>\n",
       "      <td>0.007885</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.014355</td>\n",
       "      <td>0.011524</td>\n",
       "      <td>0.020421</td>\n",
       "      <td>0.039224</td>\n",
       "      <td>0.033765</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.059644</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD1361</th>\n",
       "      <td>0.007784</td>\n",
       "      <td>0.015450</td>\n",
       "      <td>0.020050</td>\n",
       "      <td>0.038802</td>\n",
       "      <td>0.034320</td>\n",
       "      <td>0.030310</td>\n",
       "      <td>0.076660</td>\n",
       "      <td>0.040807</td>\n",
       "      <td>0.030310</td>\n",
       "      <td>0.019342</td>\n",
       "      <td>...</td>\n",
       "      <td>0.014506</td>\n",
       "      <td>0.013445</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.009317</td>\n",
       "      <td>0.012737</td>\n",
       "      <td>0.025593</td>\n",
       "      <td>0.036207</td>\n",
       "      <td>0.031608</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.068522</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD1362</th>\n",
       "      <td>0.010617</td>\n",
       "      <td>0.052726</td>\n",
       "      <td>0.009358</td>\n",
       "      <td>0.047688</td>\n",
       "      <td>0.046068</td>\n",
       "      <td>0.036890</td>\n",
       "      <td>0.064963</td>\n",
       "      <td>0.024474</td>\n",
       "      <td>0.027893</td>\n",
       "      <td>0.018355</td>\n",
       "      <td>...</td>\n",
       "      <td>0.019975</td>\n",
       "      <td>0.014216</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.011337</td>\n",
       "      <td>0.010257</td>\n",
       "      <td>0.025553</td>\n",
       "      <td>0.043549</td>\n",
       "      <td>0.028613</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.042829</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD1364</th>\n",
       "      <td>0.011400</td>\n",
       "      <td>0.024066</td>\n",
       "      <td>0.019474</td>\n",
       "      <td>0.031032</td>\n",
       "      <td>0.054782</td>\n",
       "      <td>0.031666</td>\n",
       "      <td>0.052407</td>\n",
       "      <td>0.021849</td>\n",
       "      <td>0.025966</td>\n",
       "      <td>0.015358</td>\n",
       "      <td>...</td>\n",
       "      <td>0.015833</td>\n",
       "      <td>0.015358</td>\n",
       "      <td>0.000158</td>\n",
       "      <td>0.011400</td>\n",
       "      <td>0.007600</td>\n",
       "      <td>0.023908</td>\n",
       "      <td>0.035624</td>\n",
       "      <td>0.025016</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.044332</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>LD1366</th>\n",
       "      <td>0.020267</td>\n",
       "      <td>0.042186</td>\n",
       "      <td>0.012010</td>\n",
       "      <td>0.043387</td>\n",
       "      <td>0.037682</td>\n",
       "      <td>0.051494</td>\n",
       "      <td>0.064705</td>\n",
       "      <td>0.026873</td>\n",
       "      <td>0.042036</td>\n",
       "      <td>0.019817</td>\n",
       "      <td>...</td>\n",
       "      <td>0.015613</td>\n",
       "      <td>0.013061</td>\n",
       "      <td>0.000150</td>\n",
       "      <td>0.016514</td>\n",
       "      <td>0.014412</td>\n",
       "      <td>0.051344</td>\n",
       "      <td>0.036931</td>\n",
       "      <td>0.016964</td>\n",
       "      <td>0.000000</td>\n",
       "      <td>0.044288</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>895 rows × 45 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "         TRAV1-1   TRAV1-2    TRAV10  TRAV12-1  TRAV12-2  TRAV12-3  TRAV13-1  \\\n",
       "patid                                                                          \n",
       "LD0001  0.010435  0.021739  0.016348  0.038261  0.064522  0.035130  0.057913   \n",
       "LD0002  0.007591  0.024957  0.012996  0.047039  0.043243  0.027947  0.057504   \n",
       "LD0003  0.016874  0.029601  0.011726  0.034034  0.081653  0.038467  0.062491   \n",
       "LD0006  0.013472  0.032537  0.013981  0.038638  0.039654  0.031266  0.055796   \n",
       "LD0007  0.008087  0.021781  0.012292  0.043994  0.037956  0.029761  0.074725   \n",
       "...          ...       ...       ...       ...       ...       ...       ...   \n",
       "LD1357  0.017994  0.052568  0.016175  0.038617  0.058027  0.033158  0.056814   \n",
       "LD1361  0.007784  0.015450  0.020050  0.038802  0.034320  0.030310  0.076660   \n",
       "LD1362  0.010617  0.052726  0.009358  0.047688  0.046068  0.036890  0.064963   \n",
       "LD1364  0.011400  0.024066  0.019474  0.031032  0.054782  0.031666  0.052407   \n",
       "LD1366  0.020267  0.042186  0.012010  0.043387  0.037682  0.051494  0.064705   \n",
       "\n",
       "        TRAV13-2  TRAV14/DV4    TRAV16  ...     TRAV5     TRAV6     TRAV7  \\\n",
       "patid                                   ...                                 \n",
       "LD0001  0.022957    0.032522  0.015652  ...  0.013739  0.013565  0.000000   \n",
       "LD0002  0.023117    0.021507  0.013686  ...  0.008741  0.008741  0.000000   \n",
       "LD0003  0.019305    0.037609  0.014729  ...  0.014729  0.008723  0.000000   \n",
       "LD0006  0.026563    0.021479  0.018683  ...  0.012201  0.016904  0.000000   \n",
       "LD0007  0.030192    0.039250  0.014880  ...  0.012939  0.013910  0.000000   \n",
       "...          ...         ...       ...  ...       ...       ...       ...   \n",
       "LD1357  0.024464    0.019814  0.013951  ...  0.014153  0.007885  0.000000   \n",
       "LD1361  0.040807    0.030310  0.019342  ...  0.014506  0.013445  0.000000   \n",
       "LD1362  0.024474    0.027893  0.018355  ...  0.019975  0.014216  0.000000   \n",
       "LD1364  0.021849    0.025966  0.015358  ...  0.015833  0.015358  0.000158   \n",
       "LD1366  0.026873    0.042036  0.019817  ...  0.015613  0.013061  0.000150   \n",
       "\n",
       "         TRAV8-1   TRAV8-2   TRAV8-3   TRAV8-4   TRAV8-6   TRAV9-1   TRAV9-2  \n",
       "patid                                                                         \n",
       "LD0001  0.018609  0.010435  0.029391  0.041913  0.027304  0.000000  0.071130  \n",
       "LD0002  0.011041  0.010351  0.024497  0.035998  0.027602  0.000115  0.060610  \n",
       "LD0003  0.010439  0.014443  0.021879  0.040040  0.026455  0.000000  0.047190  \n",
       "LD0006  0.015252  0.012583  0.024657  0.037367  0.027199  0.000127  0.052364  \n",
       "LD0007  0.012292  0.011861  0.026310  0.038818  0.019625  0.000000  0.074294  \n",
       "...          ...       ...       ...       ...       ...       ...       ...  \n",
       "LD1357  0.014355  0.011524  0.020421  0.039224  0.033765  0.000000  0.059644  \n",
       "LD1361  0.009317  0.012737  0.025593  0.036207  0.031608  0.000000  0.068522  \n",
       "LD1362  0.011337  0.010257  0.025553  0.043549  0.028613  0.000000  0.042829  \n",
       "LD1364  0.011400  0.007600  0.023908  0.035624  0.025016  0.000000  0.044332  \n",
       "LD1366  0.016514  0.014412  0.051344  0.036931  0.016964  0.000000  0.044288  \n",
       "\n",
       "[895 rows x 45 columns]"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "family_usage_table"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "def do_ols(pat_info_df, usage_table, ind_var, y_cols, cov_cols):\n",
    "    # ind_var: independent variable for model, e.g. 'DQB1*03:01'\n",
    "    all_df = pd.merge(pat_info_df, usage_table, right_index=True, left_index=True, how=\"inner\")\n",
    "#     if binary: # code variable as 0 1\n",
    "#         ind_var_col_name = \"{}_binary\".format(ind_var)\n",
    "#         all_df[ind_var_col_name] = (all_df[ind_var] > 0).astype(\"float\")\n",
    "#     else: # code variable as 0 1 2\n",
    "    ind_var_col_name = ind_var\n",
    "    #all_df[ind_var_col_name] = all_df[ind_var_col_name].astype(\"float\")\n",
    "    \n",
    "    #print(\"doing {} tests. out of {} rows, {} are {} > 0\".format(len(y_cols), len(all_df), len(all_df[all_df[ind_var]>0]), ind_var))\n",
    "    x_columns = [ind_var_col_name] + cov_cols\n",
    "    res_df = pd.DataFrame()\n",
    "    for y_column in y_cols:\n",
    "        X = all_df[x_columns]\n",
    "        y = all_df[y_column]\n",
    "        X = sm.add_constant(X) # constant is always added\n",
    "        mod = sm.OLS(y, X)\n",
    "        res = mod.fit()\n",
    "        #print(res.summary())\n",
    "        tcr_df = pd.DataFrame()\n",
    "        tcr_df['coef'] = res.params\n",
    "        tcr_df['se'] = res.bse\n",
    "        tcr_df['pvalue'] = res.pvalues\n",
    "        tcr_df['tvalue'] = res.tvalues\n",
    "        res_df[y_column] = tcr_df.loc[ind_var_col_name]\n",
    "    res_df = res_df.transpose().reset_index().rename(columns={\"index\":\"family\"})\n",
    "    # FDR correction\n",
    "    res_df['pvalue_fdr'] = fdrcorrection(res_df['pvalue'])[1]\n",
    "    return res_df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "def do_stepwise_ols(HLA_table, counts_table, ind_var, want_families={}, fdr=True):\n",
    "    \n",
    "    def is_continue(df):\n",
    "        if len(want_families) > 0: # need results for all families in want_families\n",
    "            print(set(df['family']).intersection(want_families), want_families)\n",
    "            if set(df['family']).intersection(want_families) == want_families:\n",
    "                return False\n",
    "        return True\n",
    "    \n",
    "    df = pd.DataFrame(columns=[\"family\", \"coef\", \"se\", \"pvalue\", \"tvalue\"])\n",
    "    drop_columns = []\n",
    "\n",
    "    while is_continue(df):\n",
    "        usage_table = counts_to_usage_table(counts_table, drop_columns=drop_columns)\n",
    "        if len(usage_table.columns) == 0: # dropped all the columns already, nothing left\n",
    "            break\n",
    "        res_df = do_ols(HLA_table, usage_table, ind_var, usage_table.columns, [])\n",
    "        top_row = res_df.iloc[0]\n",
    "        if fdr:\n",
    "            pvalue = top_row['pvalue_fdr']\n",
    "        else:\n",
    "            pvalue = top_row['pvalue']\n",
    "        if len(want_families) == 0:\n",
    "            if pvalue > 0.1: # get all rows until the top row is no longer significant\n",
    "                break\n",
    "        print(\"{}, coef: {}, pvalue: {}\".format(top_row['family'], top_row['coef'], pvalue))\n",
    "        drop_columns.append(top_row['family'])\n",
    "        df = df.append(top_row)\n",
    "    df['coef'] = df['coef']*100 # want the coef to be in percentages\n",
    "    return df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# #usage_df['usage_ratio'] = usage_df['usage_ratio']*100 # only do once!\n",
    "# cov_df = pd.read_csv(\"DGN_covariates_df.csv\").set_index(\"patid\")\n",
    "# cov_df[\"lane\"] = cov_df[\"fcid\"] + \"_\" + cov_df[\"lane\"].astype(\"str\") # covariates\n",
    "# cov_df = pd.get_dummies(cov_df, columns=['fcid', 'lane'])\n",
    "# pat_info_df = pd.merge(cov_df, HLA_table, right_index=True, left_index=True, how=\"inner\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# step_res_df_301_all = do_stepwise_ols(HLA_table, family_counts_table, ind_var=\"DQB1*03:01\", fdr=True)\n",
    "\n",
    "# step_res_df_301_all"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# tcr_columns = list(usage_table.columns)\n",
    "# do_ols(HLA_table, usage_table, ind_var=\"DQB1*03:01\", y_cols=tcr_columns, cov_cols=[])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Families OLS"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [],
   "source": [
    "for family in [\"TRAV\", \"TRBV\", \"TRAJ\", \"TRBJ\"]:\n",
    "    family_counts_table = pd.read_csv(\"DGN_family_counts_table.csv\")\n",
    "    family_counts_table = family_counts_table.set_index(\"patid\")\n",
    "    family_counts_table = family_counts_table[[x for x in family_counts_table.columns if x.startswith(family)]]\n",
    "    family_usage_table = counts_to_usage_table(family_counts_table)\n",
    "    family_res_df = do_ols(HLA_table, family_usage_table, \"DQB1*03:01\", family_usage_table.columns, [])\n",
    "    \n",
    "    # separately compute DQB0301_pos_mean_usage and DQB0301_neg_mean_usage and join in to results\n",
    "    all_df = pd.merge(HLA_table, family_usage_table, right_index=True, left_index=True, how=\"inner\")\n",
    "    num_DQB0301_pos = len(all_df[all_df[\"DQB1*03:01\"]==1.0])\n",
    "    num_DQB0301_neg = len(all_df[all_df[\"DQB1*03:01\"]!=1.0])\n",
    "    DQB0301_pos_mean_usage = all_df[all_df[\"DQB1*03:01\"]==1.0][[x for x in family_counts_table.columns if x.startswith(family)]].mean().to_frame(name=f\"DQB0301_pos_mean_usage\").reset_index().rename(columns={'index': 'family'})\n",
    "    DQB0301_neg_mean_usage = all_df[all_df[\"DQB1*03:01\"]!=1.0][[x for x in family_counts_table.columns if x.startswith(family)]].mean().to_frame(name=f\"DQB0301_neg_mean_usage\").reset_index().rename(columns={'index': 'family'})\n",
    "    mean_usage_df = DQB0301_pos_mean_usage.merge(DQB0301_neg_mean_usage, how=\"inner\", on=\"family\")\n",
    "    mean_usage_df[\"DQB0301_pos_to_neg_mean_usage_ratio\"] = mean_usage_df[\"DQB0301_pos_mean_usage\"] / mean_usage_df[\"DQB0301_neg_mean_usage\"]\n",
    "    merged_df = mean_usage_df.merge(family_res_df, how=\"inner\", on=\"family\")\n",
    "    merged_df = merged_df.sort_values(by=['pvalue_fdr'])\n",
    "    merged_df = merged_df.rename(columns={\"DQB0301_pos_mean_usage\": f\"DQB0301_pos_mean_usage (n={num_DQB0301_pos})\", \"DQB0301_neg_mean_usage\": f\"DQB0301_neg_mean_usage (n={num_DQB0301_neg})\"})\n",
    "    merged_df.to_csv(f\"FINAL_DGN_{family}_FAMILY_OLS.csv\", sep=\"\\t\", index=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Families Stepwise OLS"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "step_family_res_df = do_stepwise_ols(HLA_table, family_counts_table, \"DQB1*03:01\")\n",
    "step_family_res_df.to_csv(\"FINAL_DGN_FAMILY_STEPWISE_OLS.csv\", sep=\"\\t\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cond_A = HLA_table['DQB1*03:01'] > 0\n",
    "cond_B = HLA_table['DQB1*03:02'] > 0\n",
    "filt_HLA_table = HLA_table[(cond_A & ~cond_B) | (cond_B & ~cond_A)]\n",
    "step_res_df_301_302 = do_stepwise_ols(filt_HLA_table, family_counts_table, \"DQB1*03:01\", fdr=True, covariates=False, binary=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "step_res_df_301_302"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cond_A = HLA_table['DQB1*03:01'] > 0\n",
    "cond_B = HLA_table['DQB1*06:01'] > 0\n",
    "filt_HLA_table = HLA_table[(cond_A & ~cond_B) | (cond_B & ~cond_A)]\n",
    "step_res_df_301_302 = do_stepwise_ols(filt_HLA_table, family_counts_table, \"DQB1*03:01\", fdr=True, covariates=False, binary=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Amino Acids OLS"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "ind_var = \"DQB1*03:01\"\n",
    "#print(\"indicator variable is {}\".format(ind_var))\n",
    "\n",
    "res_df_lst = []\n",
    "\n",
    "## CDR1\n",
    "TRAV_active_positions = [0, 1, 2, 3, 9, 10, 11]\n",
    "TRBV_active_positions = [0, 1, 2, 9, 10, 11]\n",
    "\n",
    "if TRV == \"TRAV\":\n",
    "    TRV_active_positions = TRAV_active_positions\n",
    "else:\n",
    "    TRV_active_positions = TRBV_active_positions\n",
    "\n",
    "#print(TRV, \"CDR1\")\n",
    "#print(\"everyone with {} is also {}\".format(ind_var, filt_col))\n",
    "for i in TRV_active_positions:\n",
    "    #print(\"--- position {} ---\".format(i+27))\n",
    "    position_AA_usage_table = counts_to_usage_table(positions_AA_count_dfs[i])\n",
    "    res_df = do_ols(HLA_table, position_AA_usage_table, ind_var, position_AA_usage_table.columns, [])\n",
    "    res_df[\"position\"] = i+27\n",
    "    res_df_lst.append(res_df)\n",
    "\n",
    "#print(\"\\n\")\n",
    "## CDR2\n",
    "TRAV_active_positions = [0, 1, 2, 3, 6, 7, 8, 9]\n",
    "TRBV_active_positions = [0, 1, 2, 3, 7, 8, 9]\n",
    "\n",
    "if TRV == \"TRAV\":\n",
    "    TRV_active_positions = TRAV_active_positions\n",
    "else:\n",
    "    TRV_active_positions = TRBV_active_positions\n",
    "\n",
    "#print(TRV, \"CDR2\")\n",
    "#print(\"everyone with {} is also {}\".format(ind_var, filt_col))\n",
    "for i in TRV_active_positions:\n",
    "    #print(\"--- position {} ---\".format(i+56))\n",
    "    # +12 to offset the CDR1 AA's in front\n",
    "    position_AA_usage_table = counts_to_usage_table(positions_AA_count_dfs[i+12])\n",
    "    res_df = do_ols(HLA_table, position_AA_usage_table, ind_var, position_AA_usage_table.columns, [])\n",
    "    res_df[\"position\"] = i+56\n",
    "    res_df_lst.append(res_df)\n",
    "    \n",
    "chi_square_AA_res_df = pd.concat(res_df_lst)\n",
    "chi_square_AA_res_df = chi_square_AA_res_df[[\"position\", \"family\", \"coef\", \"pvalue\"]]\n",
    "chi_square_AA_res_df['pvalue_fdr'] = fdrcorrection(chi_square_AA_res_df['pvalue'])[1]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "chi_square_AA_res_df[chi_square_AA_res_df[\"pvalue_fdr\"] < 0.05].sort_values(by=['position', 'pvalue_fdr']).to_csv(f\"FINAL_DGN_{TRV}_AMINOACID_OLS.csv\", sep=\"\\t\", index=None)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Amino Acids Stepwise OLS"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#print(\"everyone is either DQA1*03:03 or DQA1*05:05\")\n",
    "ind_var = \"DQB1*03:01\" #\"DQA1*05:05\" \n",
    "print(\"indicator variable is {}\".format(ind_var))\n",
    "# for filt_col in ['DQA1*05:05', 'DQA1*03:03']:\n",
    "print(\"\\n\")\n",
    "#HLA_table_filt = HLA_table[(HLA_table['DQA1*03:03'] == 1) | (HLA_table['DQA1*05:05'] == 1)]\n",
    "# HLA_table_filt = HLA_table[~((HLA_table['DQB1*03:01'] == 1) & (HLA_table[filt_col] != 1))]\n",
    "# HLA_table_filt = HLA_table[HLA_table[filt_col] == 1.0] # HLA_table[HLA_table['DQA1*03:03'] == 1.0]\n",
    "HLA_table_filt = HLA_table\n",
    "\n",
    "## CDR1\n",
    "TRAV_active_positions = [0, 1, 2, 3, 9, 10, 11]\n",
    "TRBV_active_positions = [0, 1, 2, 9, 10, 11]\n",
    "\n",
    "if TRV == \"TRAV\":\n",
    "    TRV_active_positions = TRAV_active_positions\n",
    "else:\n",
    "    TRV_active_positions = TRBV_active_positions\n",
    "\n",
    "print(TRV, \"CDR1\")\n",
    "#print(\"everyone with {} is also {}\".format(ind_var, filt_col))\n",
    "for i in TRV_active_positions:\n",
    "    print(\"--- position {} ---\".format(i+27))\n",
    "    do_stepwise_ols(HLA_table_filt, positions_AA_count_dfs[i], ind_var=ind_var, fdr=True)\n",
    "\n",
    "print(\"\\n\")\n",
    "## CDR2\n",
    "TRAV_active_positions = [0, 1, 2, 3, 6, 7, 8, 9]\n",
    "TRBV_active_positions = [0, 1, 2, 3, 7, 8, 9]\n",
    "\n",
    "if TRV == \"TRAV\":\n",
    "    TRV_active_positions = TRAV_active_positions\n",
    "else:\n",
    "    TRV_active_positions = TRBV_active_positions\n",
    "\n",
    "print(TRV, \"CDR2\")\n",
    "#print(\"everyone with {} is also {}\".format(ind_var, filt_col))\n",
    "for i in TRV_active_positions:\n",
    "    print(\"--- position {} ---\".format(i+56))\n",
    "    # +12 to offset the CDR1 AA's in front\n",
    "    do_stepwise_ols(HLA_table_filt, positions_AA_count_dfs[i+12], ind_var=ind_var, fdr=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "positions_AA_count_dfs[0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python3.7 (fsl)",
   "language": "python",
   "name": "fsl"
  },
  "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.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}