{
"cells": [
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
"# Inline graphing\n",
"%matplotlib inline"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"# Import necessary packages\n",
"import os\n",
"import pandas as pd\n",
"import numpy as np\n",
"import seaborn as sns\n",
"import matplotlib.pyplot as plt\n",
"from matplotlib_venn import venn3, venn2\n",
"from scipy import stats\n",
"\n",
"pd.set_option('mode.chained_assignment', None)"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [],
"source": [
"# Set working directory to wild-type 1 folder\n",
"os.chdir('Path/To/Data')"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [],
"source": [
"# Reading in full data for first wild-type experiment\n",
"ua1 = pd.read_csv('UnA/rowsum/20190516_UnA_SVCidentified_postiter_threshold.csv', index_col=0)\n",
"ua1.drop(columns=['3K', '5.4K', '12.2K', '24K', '78.4K', '110K', '195.5K'], inplace=True)\n",
"ua1.columns = ['Gene_ua1', 'ProteinInfo_ua1', 'Compartment_ua1', 'Probability_ua1']\n",
"\n",
"ub1 = pd.read_csv('UnB/rowsum/20190516_UnB_SVCidentified_postiter_threshold.csv', index_col=0)\n",
"ub1.drop(columns=['3K', '5.4K', '12.2K', '24K', '78.4K', '110K', '195.5K'], inplace=True)\n",
"ub1.columns = ['Gene_ub1', 'ProteinInfo_ub1', 'Compartment_ub1', 'Probability_ub1']\n",
"\n",
"uc1 = pd.read_csv('UnC/rowsum/20190516_UnC_SVCidentified_postiter_threshold.csv', index_col=0)\n",
"uc1.drop(columns=['3K', '5.4K', '12.2K', '24K', '78.4K', '110K', '195.5K'], inplace=True)\n",
"uc1.columns = ['Gene_uc1', 'ProteinInfo_uc1', 'Compartment_uc1', 'Probability_uc1']\n",
"\n",
"ia1 = pd.read_csv('IndA/rowsum/20190516_IndA_SVCidentified_postiter_threshold.csv', index_col=0)\n",
"ia1.drop(columns=['3K', '5.4K', '12.2K', '24K', '78.4K', '110K', '195.5K'], inplace=True)\n",
"ia1.columns = ['Gene_ia1', 'ProteinInfo_ia1', 'Compartment_ia1', 'Probability_ia1']\n",
"\n",
"ib1 = pd.read_csv('IndB/rowsum/20190516_IndB_SVCidentified_postiter_threshold.csv', index_col=0)\n",
"ib1.drop(columns=['3K', '5.4K', '12.2K', '24K', '78.4K', '110K', '195.5K'], inplace=True)\n",
"ib1.columns = ['Gene_ib1', 'ProteinInfo_ib1', 'Compartment_ib1', 'Probability_ib1']\n",
"\n",
"ic1 = pd.read_csv('IndC/rowsum/20190516_IndC_SVCidentified_postiter_threshold.csv', index_col=0)\n",
"ic1.drop(columns=['3K', '5.4K', '12.2K', '24K', '78.4K', '110K', '195.5K'], inplace=True)\n",
"ic1.columns = ['Gene_ic1', 'ProteinInfo_ic1', 'Compartment_ic1', 'Probability_ic1']"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [],
"source": [
"# Set working directory to wild-type 2 folder\n",
"os.chdir('Path/To/Data')"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [],
"source": [
"# Reading in full data for second wild-type experiment\n",
"ua2 = pd.read_csv('UnA/20190518_UnA_SVCidentified_postiter_threshold.csv', index_col=0)\n",
"ua2.drop(columns=['3K', '5.4K', '12.2K', '24K', '78.4K', '110K', '195.5K'], inplace=True)\n",
"ua2.columns = ['Gene_ua2', 'ProteinInfo_ua2', 'Compartment_ua2', 'Probability_ua2']\n",
"\n",
"ub2 = pd.read_csv('UnB/20190518_UnB_SVCidentified_postiter_threshold.csv', index_col=0)\n",
"ub2.drop(columns=['3K', '5.4K', '12.2K', '24K', '78.4K', '110K', '195.5K'], inplace=True)\n",
"ub2.columns = ['Gene_ub2', 'ProteinInfo_ub2', 'Compartment_ub2', 'Probability_ub2']\n",
"\n",
"uc2 = pd.read_csv('UnC/20190518_UnC_SVCidentified_postiter_threshold.csv', index_col=0)\n",
"uc2.drop(columns=['3K', '5.4K', '12.2K', '24K', '78.4K', '110K', '195.5K'], inplace=True)\n",
"uc2.columns = ['Gene_uc2', 'ProteinInfo_uc2', 'Compartment_uc2', 'Probability_uc2']\n",
"\n",
"ia2 = pd.read_csv('IndA/20190518_IndA_SVCidentified_postiter_threshold.csv', index_col=0)\n",
"ia2.drop(columns=['3K', '5.4K', '12.2K', '24K', '78.4K', '110K', '195.5K'], inplace=True)\n",
"ia2.columns = ['Gene_ia2', 'ProteinInfo_ia2', 'Compartment_ia2', 'Probability_ia2']\n",
"\n",
"ib2 = pd.read_csv('IndB/20190518_IndB_SVCidentified_postiter_threshold.csv', index_col=0)\n",
"ib2.drop(columns=['3K', '5.4K', '12.2K', '24K', '78.4K', '110K', '195.5K'], inplace=True)\n",
"ib2.columns = ['Gene_ib2', 'ProteinInfo_ib2', 'Compartment_ib2', 'Probability_ib2']\n",
"\n",
"ic2 = pd.read_csv('IndC/20190518_IndC_SVCidentified_postiter_threshold.csv', index_col=0)\n",
"ic2.drop(columns=['3K', '5.4K', '12.2K', '24K', '78.4K', '110K', '195.5K'], inplace=True)\n",
"ic2.columns = ['Gene_ic2', 'ProteinInfo_ic2', 'Compartment_ic2', 'Probability_ic2']"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [],
"source": [
"# Set working directory to folder comparing biological replicates\n",
"os.chdir('Path/To/Data')"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [],
"source": [
"# Defining commonly detected protein dataframes for uninduced\n",
"unlist = [ua1, ub1, uc1, ua2, ub2, uc2]\n",
"\n",
"un = pd.concat(unlist, axis=1, join='outer', sort='True')\n",
"\n",
"for item in un.index:\n",
" if un.at[item, 'ProteinInfo_ua1'] is np.nan:\n",
" if un.at[item, 'ProteinInfo_ub1'] is not np.nan:\n",
" un.at[item, 'ProteinInfo_ua1'] = un.at[item, 'ProteinInfo_ub1']\n",
" elif un.at[item, 'ProteinInfo_uc1'] is not np.nan:\n",
" un.at[item, 'ProteinInfo_ua1'] = un.at[item, 'ProteinInfo_uc1']\n",
" elif un.at[item, 'ProteinInfo_ua2'] is not np.nan:\n",
" un.at[item, 'ProteinInfo_ua1'] = un.at[item, 'ProteinInfo_ua2']\n",
" elif un.at[item, 'ProteinInfo_ub2'] is not np.nan:\n",
" un.at[item, 'ProteinInfo_ua1'] = un.at[item, 'ProteinInfo_ub2']\n",
" elif un.at[item, 'ProteinInfo_uc2'] is not np.nan:\n",
" un.at[item, 'ProteinInfo_ua1'] = un.at[item, 'ProteinInfo_uc2']\n",
" \n",
"for item in un.index:\n",
" if un.at[item, 'Gene_ua1'] is np.nan:\n",
" if un.at[item, 'Gene_ub1'] is not np.nan:\n",
" un.at[item, 'Gene_ua1'] = un.at[item, 'Gene_ub1']\n",
" elif un.at[item, 'Gene_uc1'] is not np.nan:\n",
" un.at[item, 'Gene_ua1'] = un.at[item, 'Gene_uc1']\n",
" elif un.at[item, 'Gene_ua2'] is not np.nan:\n",
" un.at[item, 'Gene_ua1'] = un.at[item, 'Gene_ua2']\n",
" elif un.at[item, 'Gene_ub2'] is not np.nan:\n",
" un.at[item, 'Gene_ua1'] = un.at[item, 'Gene_ub2']\n",
" elif un.at[item, 'Gene_uc2'] is not np.nan:\n",
" un.at[item, 'Gene_ua1'] = un.at[item, 'Gene_uc2']\n",
"\n",
"un.drop(columns=['ProteinInfo_ub1', 'ProteinInfo_uc1', 'ProteinInfo_ua2', 'ProteinInfo_ub2', 'ProteinInfo_uc2', \n",
" 'Gene_ub1', 'Gene_uc1', 'Gene_ua2', 'Gene_ub2', 'Gene_uc2'], inplace=True)\n",
"un.rename({'Gene_ua1':'Gene', 'ProteinInfo_ua1':'ProteinInfo'}, axis='columns', inplace=True)\n",
"un.fillna('ND', inplace=True)\n",
"ndset = ['ND']\n",
"un['NDcount'] = (un.isin(ndset).sum(1))/2"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {},
"outputs": [],
"source": [
"# Defining commonly detected protein dataframes for induced\n",
"indlist = [ia1, ib1, ic1, ia2, ib2, ic2]\n",
"\n",
"ind = pd.concat(indlist, axis=1, join='outer', sort='True')\n",
"\n",
"for item in ind.index:\n",
" if ind.at[item, 'ProteinInfo_ia1'] is np.nan:\n",
" if ind.at[item, 'ProteinInfo_ib1'] is not np.nan:\n",
" ind.at[item, 'ProteinInfo_ia1'] = ind.at[item, 'ProteinInfo_ib1']\n",
" elif ind.at[item, 'ProteinInfo_ic1'] is not np.nan:\n",
" ind.at[item, 'ProteinInfo_ia1'] = ind.at[item, 'ProteinInfo_ic1']\n",
" elif ind.at[item, 'ProteinInfo_ia2'] is not np.nan:\n",
" ind.at[item, 'ProteinInfo_ia1'] = ind.at[item, 'ProteinInfo_ia2']\n",
" elif ind.at[item, 'ProteinInfo_ib2'] is not np.nan:\n",
" ind.at[item, 'ProteinInfo_ia1'] = ind.at[item, 'ProteinInfo_ib2']\n",
" elif ind.at[item, 'ProteinInfo_ic2'] is not np.nan:\n",
" ind.at[item, 'ProteinInfo_ia1'] = ind.at[item, 'ProteinInfo_ic2']\n",
" \n",
"for item in ind.index:\n",
" if ind.at[item, 'Gene_ia1'] is np.nan:\n",
" if ind.at[item, 'Gene_ib1'] is not np.nan:\n",
" ind.at[item, 'Gene_ia1'] = ind.at[item, 'Gene_ib1']\n",
" elif ind.at[item, 'Gene_ic1'] is not np.nan:\n",
" ind.at[item, 'Gene_ia1'] = ind.at[item, 'Gene_ic1']\n",
" elif ind.at[item, 'Gene_ia2'] is not np.nan:\n",
" ind.at[item, 'Gene_ia1'] = ind.at[item, 'Gene_ia2']\n",
" elif ind.at[item, 'Gene_ib2'] is not np.nan:\n",
" ind.at[item, 'Gene_ia1'] = ind.at[item, 'Gene_ib2']\n",
" elif ind.at[item, 'Gene_ic2'] is not np.nan:\n",
" ind.at[item, 'Gene_ia1'] = ind.at[item, 'Gene_ic2']\n",
"\n",
"ind.drop(columns=['ProteinInfo_ib1', 'ProteinInfo_ic1', 'ProteinInfo_ia2', 'ProteinInfo_ib2', 'ProteinInfo_ic2', \n",
" 'Gene_ib1', 'Gene_ic1', 'Gene_ia2', 'Gene_ib2', 'Gene_ic2'], inplace=True)\n",
"ind.rename({'Gene_ia1':'Gene', 'ProteinInfo_ia1':'ProteinInfo'}, axis='columns', inplace=True)\n",
"ind.fillna('ND', inplace=True)\n",
"ndset = ['ND']\n",
"ind['NDcount'] = (ind.isin(ndset).sum(1))/2"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {},
"outputs": [],
"source": [
"# Creating dataframes with only those proteins classified in all 6 replicates for a given condition (excludes proteins that were markers for\n",
"# both WT biological replicates)\n",
"un_common = un[un['NDcount'] == 0.0]\n",
"un_common = un_common[((un_common['Probability_ua1'] > 0) & (un_common['Probability_ua2'] > 0))]\n",
"\n",
"ind_common = ind[ind['NDcount'] == 0.0]\n",
"ind_common = ind_common[((ind_common['Probability_ia1'] > 0) & (ind_common['Probability_ia2'] > 0))]"
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {},
"outputs": [],
"source": [
"# Generating dataframes for finding row modes in compartment columns\n",
"un_temp = un_common.drop(columns=['Probability_ua1', 'Probability_ub1', 'Probability_uc1', \n",
" 'Probability_ua2', 'Probability_ub2', 'Probability_uc2'])\n",
"un_temp = un_temp.loc[:, 'Compartment_ua1':'Compartment_uc2']\n",
"\n",
"ind_temp = ind_common.drop(columns=['Probability_ia1', 'Probability_ib1', 'Probability_ic1', \n",
" 'Probability_ia2', 'Probability_ib2', 'Probability_ic2'])\n",
"ind_temp = ind_temp.loc[:, 'Compartment_ia1':'Compartment_ic2']\n",
"\n",
"# Taking row mode and adding this data back onto the dataframes with proteins classified in all 6 replicates for a given condition\n",
"a = un_temp.mode(axis=1)\n",
"a.columns = ['Mode1', 'Mode2', 'Mode3']\n",
"un_common = pd.concat([un_common, a], axis=1, join='inner')\n",
"\n",
"b = ind_temp.mode(axis=1)\n",
"b.columns = ['Mode1', 'Mode2', 'Mode3']\n",
"ind_common = pd.concat([ind_common, b], axis=1, join='inner')"
]
},
{
"cell_type": "code",
"execution_count": 12,
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"C:\\Users\\Aaron Oom\\Anaconda3\\lib\\site-packages\\scipy\\stats\\stats.py:245: RuntimeWarning: The input array could not be properly checked for nan values. nan values will be ignored.\n",
" \"values. nan values will be ignored.\", RuntimeWarning)\n"
]
}
],
"source": [
"# Counting instances of row mode for each protein\n",
"x, y = stats.mode(un_temp, axis=1)\n",
"un_common['ModeCount'] = y\n",
"\n",
"x, y = stats.mode(ind_temp, axis=1)\n",
"ind_common['ModeCount'] = y"
]
},
{
"cell_type": "code",
"execution_count": 13,
"metadata": {},
"outputs": [],
"source": [
"# Generating dataframes for each condition with at least 4, 5, or 6 replicates classified the same\n",
"un6 = un_common[un_common.ModeCount >= 6]\n",
"un6.drop(columns=['NDcount', 'Mode2', 'Mode3'], inplace=True)\n",
"un6.rename({'Mode1':'Mode_un'}, inplace=True, axis='columns')\n",
"\n",
"un5 = un_common[un_common.ModeCount >= 5]\n",
"un5.drop(columns=['NDcount', 'Mode2', 'Mode3'], inplace=True)\n",
"un5.rename({'Mode1':'Mode_un'}, inplace=True, axis='columns')\n",
"\n",
"un4 = un_common[un_common.ModeCount >= 4]\n",
"un4.drop(columns=['NDcount', 'Mode2', 'Mode3'], inplace=True)\n",
"un4.rename({'Mode1':'Mode_un'}, inplace=True, axis='columns')\n",
"\n",
"ind6 = ind_common[ind_common.ModeCount >= 6]\n",
"ind6.drop(columns=['NDcount', 'Mode2', 'Mode3', 'Gene', 'ProteinInfo'], inplace=True)\n",
"ind6.rename({'Mode1':'Mode_ind'}, inplace=True, axis='columns')\n",
"\n",
"ind5 = ind_common[ind_common.ModeCount >= 5]\n",
"ind5.drop(columns=['NDcount', 'Mode2', 'Mode3', 'Gene', 'ProteinInfo'], inplace=True)\n",
"ind5.rename({'Mode1':'Mode_ind'}, inplace=True, axis='columns')\n",
"\n",
"ind4 = ind_common[ind_common.ModeCount >= 4]\n",
"ind4.drop(columns=['NDcount', 'Mode2', 'Mode3', 'Gene', 'ProteinInfo'], inplace=True)\n",
"ind4.rename({'Mode1':'Mode_ind'}, inplace=True, axis='columns')"
]
},
{
"cell_type": "code",
"execution_count": 14,
"metadata": {},
"outputs": [],
"source": [
"# Concatenating uninduced and induced for each cutoff\n",
"common6 = pd.concat([un6, ind6], join='inner', axis=1, sort=True)\n",
"common5 = pd.concat([un5, ind5], join='inner', axis=1, sort=True)\n",
"common4 = pd.concat([un4, ind4], join='inner', axis=1, sort=True)"
]
},
{
"cell_type": "code",
"execution_count": 15,
"metadata": {},
"outputs": [],
"source": [
"# Subsetting out proteins that have moved\n",
"common6_moved = common6[common6['Mode_un'] != common6['Mode_ind']]\n",
"common5_moved = common5[common5['Mode_un'] != common5['Mode_ind']]\n",
"common4_moved = common4[common4['Mode_un'] != common4['Mode_ind']]"
]
},
{
"cell_type": "code",
"execution_count": 16,
"metadata": {},
"outputs": [],
"source": [
"# Saving common5_moved and common4_moved\n",
"common5_moved.to_csv('20190528_RelocalizedProteins_5of6repl.csv', sep=',')\n",
"common5_ID = common5_moved.loc[:, 'Gene':'ProteinInfo']\n",
"common5_ID.reset_index(inplace=True)\n",
"common5_ID.drop(['Gene', 'ProteinInfo'], axis=1, inplace=True)\n",
"common5_ID.to_csv('20190528_RelocalizedProteins_5of6repl_IDs.csv', sep=',', index=False)\n",
"\n",
"common4_moved.to_csv('20190528_RelocalizedProteins_4of6repl.csv', sep=',')\n",
"common4_ID = common4_moved.loc[:, 'Gene':'ProteinInfo']\n",
"common4_ID.reset_index(inplace=True)\n",
"common4_ID.drop(['Gene', 'ProteinInfo'], axis=1, inplace=True)\n",
"common4_ID.to_csv('20190528_RelocalizedProteins_4of6repl_IDs.csv', sep=',', index=False)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### STRING Analyses for Moved Proteins\n",
"\n",
"[Proteins identical in 4+ replicates/condition](https://version-11-0.string-db.org/cgi/network.pl?networkId=0pi4Kj1XKAjO)<br>\n",
"[Proteins identical in 5+ replicates/condition](https://version-11-0.string-db.org/cgi/network.pl?networkId=HBuAuan2xdNC)<br>"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Centroid Distance Analysis"
]
},
{
"cell_type": "code",
"execution_count": 57,
"metadata": {},
"outputs": [],
"source": [
"# Import packages\n",
"import os\n",
"import math\n",
"import pandas as pd\n",
"import numpy as np\n",
"import seaborn as sns\n",
"import matplotlib.pyplot as plt\n",
"from sklearn.covariance import MinCovDet\n",
"from scipy.stats import percentileofscore\n",
"from scipy.stats import rankdata\n",
"from scipy.stats import t\n",
"from scipy.spatial.distance import euclidean"
]
},
{
"cell_type": "code",
"execution_count": 58,
"metadata": {},
"outputs": [],
"source": [
"# Set working directory to wild-type 1 folder\n",
"os.chdir('C:/Users/Aaron Oom/Documents/UCSD/Guatelli Lab/Cell Fractionation/20190308 Jurkat TREHIV WT fract MS/')"
]
},
{
"cell_type": "code",
"execution_count": 59,
"metadata": {},
"outputs": [],
"source": [
"# Reading in full data for first wild-type experiment\n",
"ua1 = pd.read_csv('UnA/rowsum/20190516_UnA_RowSum.csv', index_col=0)\n",
"ub1 = pd.read_csv('UnB/rowsum/20190516_UnB_RowSum.csv', index_col=0)\n",
"uc1 = pd.read_csv('UnC/rowsum/20190516_UnC_RowSum.csv', index_col=0)\n",
"\n",
"ia1 = pd.read_csv('IndA/rowsum/20190516_IndA_RowSum.csv', index_col=0)\n",
"ib1 = pd.read_csv('IndB/rowsum/20190516_IndB_RowSum.csv', index_col=0)\n",
"ic1 = pd.read_csv('IndC/rowsum/20190516_IndC_RowSum.csv', index_col=0)\n",
"\n",
"wt1_moved = pd.read_csv('20190516_EuclideanDistance_FDR0.02.csv', index_col=0)"
]
},
{
"cell_type": "code",
"execution_count": 60,
"metadata": {},
"outputs": [],
"source": [
"# Labeling replicate columns for ease of identification in joined dataframe\n",
"ua1.columns = ['Gene_ua1', 'Protein information_ua1', '3K_ua1', '5.4K_ua1', '12.2K_ua1', \n",
" '24K_ua1', '78.4K_ua1', '110K_ua1', '195.5K_ua1']\n",
"ub1.columns = ['Gene_ub1', 'Protein information_ub1', '3K_ub1', '5.4K_ub1', '12.2K_ub1', \n",
" '24K_ub1', '78.4K_ub1', '110K_ub1', '195.5K_ub1']\n",
"uc1.columns = ['Gene_uc1', 'Protein information_uc1', '3K_uc1', '5.4K_uc1', '12.2K_uc1', \n",
" '24K_uc1', '78.4K_uc1', '110K_uc1', '195.5K_uc1']\n",
"ia1.columns = ['Gene_ia1', 'Protein information_ia1', '3K_ia1', '5.4K_ia1', '12.2K_ia1', \n",
" '24K_ia1', '78.4K_ia1', '110K_ia1', '195.5K_ia1']\n",
"ib1.columns = ['Gene_ib1', 'Protein information_ib1', '3K_ib1', '5.4K_ib1', '12.2K_ib1', \n",
" '24K_ib1', '78.4K_ib1', '110K_ib1', '195.5K_ib1']\n",
"ic1.columns = ['Gene_ic1', 'Protein information_ic1', '3K_ic1', '5.4K_ic1', '12.2K_ic1', \n",
" '24K_ic1', '78.4K_ic1', '110K_ic1', '195.5K_ic1']"
]
},
{
"cell_type": "code",
"execution_count": 61,
"metadata": {},
"outputs": [],
"source": [
"# Set working directory to wild-type 2 folder\n",
"os.chdir('C:/Users/Aaron Oom/Documents/UCSD/Guatelli Lab/Cell Fractionation/20190517 Jurkat TREHIV WT and dNef fract MS/')"
]
},
{
"cell_type": "code",
"execution_count": 62,
"metadata": {},
"outputs": [],
"source": [
"# Reading in full data for second wild-type experiment\n",
"ua2 = pd.read_csv('UnA/20190517_UnA_RowSum.csv', index_col=0)\n",
"ub2 = pd.read_csv('UnB/20190517_UnB_RowSum.csv', index_col=0)\n",
"uc2 = pd.read_csv('UnC/20190517_UnC_RowSum.csv', index_col=0)\n",
"\n",
"ia2 = pd.read_csv('IndA/20190517_IndA_RowSum.csv', index_col=0)\n",
"ib2 = pd.read_csv('IndB/20190517_IndB_RowSum.csv', index_col=0)\n",
"ic2 = pd.read_csv('IndC/20190517_IndC_RowSum.csv', index_col=0)\n",
"\n",
"wt2_moved = pd.read_csv('20190518_EuclideanDistance_FDR0.015.csv', index_col=0)"
]
},
{
"cell_type": "code",
"execution_count": 63,
"metadata": {},
"outputs": [],
"source": [
"# Labeling replicate columns for ease of identification in joined dataframe\n",
"ua2.columns = ['Gene_ua2', 'Protein information_ua2', '3K_ua2', '5.4K_ua2', '12.2K_ua2', \n",
" '24K_ua2', '78.4K_ua2', '110K_ua2', '195.5K_ua2']\n",
"ub2.columns = ['Gene_ub2', 'Protein information_ub2', '3K_ub2', '5.4K_ub2', '12.2K_ub2', \n",
" '24K_ub2', '78.4K_ub2', '110K_ub2', '195.5K_ub2']\n",
"uc2.columns = ['Gene_uc2', 'Protein information_uc2', '3K_uc2', '5.4K_uc2', '12.2K_uc2', \n",
" '24K_uc2', '78.4K_uc2', '110K_uc2', '195.5K_uc2']\n",
"ia2.columns = ['Gene_ia2', 'Protein information_ia2', '3K_ia2', '5.4K_ia2', '12.2K_ia2', \n",
" '24K_ia2', '78.4K_ia2', '110K_ia2', '195.5K_ia2']\n",
"ib2.columns = ['Gene_ib2', 'Protein information_ib2', '3K_ib2', '5.4K_ib2', '12.2K_ib2', \n",
" '24K_ib2', '78.4K_ib2', '110K_ib2', '195.5K_ib2']\n",
"ic2.columns = ['Gene_ic2', 'Protein information_ic2', '3K_ic2', '5.4K_ic2', '12.2K_ic2', \n",
" '24K_ic2', '78.4K_ic2', '110K_ic2', '195.5K_ic2']"
]
},
{
"cell_type": "code",
"execution_count": 64,
"metadata": {},
"outputs": [],
"source": [
"# Set working directory to folder comparing biological replicates\n",
"os.chdir('C:/Users/Aaron Oom/Documents/UCSD/Guatelli Lab/Cell Fractionation/20190517 Jurkat TREHIV WT and dNef fract MS/ComparingBiolRepl')"
]
},
{
"cell_type": "code",
"execution_count": 65,
"metadata": {},
"outputs": [],
"source": [
"# Pull out proteins common across all six replicates\n",
"dflist = [ua1, ub1, uc1, ia1, ib1, ic1, ua2, ub2, uc2, ia2, ib2, ic2]\n",
"common = pd.concat(dflist, axis=1, join='inner')"
]
},
{
"cell_type": "code",
"execution_count": 66,
"metadata": {},
"outputs": [],
"source": [
"# Redefine replicates with only common proteins\n",
"ua1 = common.loc[:, 'Gene_ua1':'195.5K_ua1']\n",
"ub1 = common.loc[:, 'Gene_ub1':'195.5K_ub1']\n",
"uc1 = common.loc[:, 'Gene_uc1':'195.5K_uc1']\n",
"ia1 = common.loc[:, 'Gene_ia1':'195.5K_ia1']\n",
"ib1 = common.loc[:, 'Gene_ib1':'195.5K_ib1']\n",
"ic1 = common.loc[:, 'Gene_ic1':'195.5K_ic1']\n",
"\n",
"ua2 = common.loc[:, 'Gene_ua2':'195.5K_ua2']\n",
"ub2 = common.loc[:, 'Gene_ub2':'195.5K_ub2']\n",
"uc2 = common.loc[:, 'Gene_uc2':'195.5K_uc2']\n",
"ia2 = common.loc[:, 'Gene_ia2':'195.5K_ia2']\n",
"ib2 = common.loc[:, 'Gene_ib2':'195.5K_ib2']\n",
"ic2 = common.loc[:, 'Gene_ic2':'195.5K_ic2']"
]
},
{
"cell_type": "code",
"execution_count": 67,
"metadata": {},
"outputs": [],
"source": [
"# Changing columns back to replicate invariant names\n",
"ua1.columns = ['Gene', 'Protein information', '3K', '5.4K', '12.2K', '24K', '78.4K', '110K', '195.5K']\n",
"ub1.columns = ['Gene', 'Protein information', '3K', '5.4K', '12.2K', '24K', '78.4K', '110K', '195.5K']\n",
"uc1.columns = ['Gene', 'Protein information', '3K', '5.4K', '12.2K', '24K', '78.4K', '110K', '195.5K']\n",
"ia1.columns = ['Gene', 'Protein information', '3K', '5.4K', '12.2K', '24K', '78.4K', '110K', '195.5K']\n",
"ib1.columns = ['Gene', 'Protein information', '3K', '5.4K', '12.2K', '24K', '78.4K', '110K', '195.5K']\n",
"ic1.columns = ['Gene', 'Protein information', '3K', '5.4K', '12.2K', '24K', '78.4K', '110K', '195.5K']\n",
"\n",
"ua2.columns = ['Gene', 'Protein information', '3K', '5.4K', '12.2K', '24K', '78.4K', '110K', '195.5K']\n",
"ub2.columns = ['Gene', 'Protein information', '3K', '5.4K', '12.2K', '24K', '78.4K', '110K', '195.5K']\n",
"uc2.columns = ['Gene', 'Protein information', '3K', '5.4K', '12.2K', '24K', '78.4K', '110K', '195.5K']\n",
"ia2.columns = ['Gene', 'Protein information', '3K', '5.4K', '12.2K', '24K', '78.4K', '110K', '195.5K']\n",
"ib2.columns = ['Gene', 'Protein information', '3K', '5.4K', '12.2K', '24K', '78.4K', '110K', '195.5K']\n",
"ic2.columns = ['Gene', 'Protein information', '3K', '5.4K', '12.2K', '24K', '78.4K', '110K', '195.5K']"
]
},
{
"cell_type": "code",
"execution_count": 68,
"metadata": {},
"outputs": [],
"source": [
"# Calculate centroid arrays for uninduced and induced\n",
"un_cent = np.zeros((len(common.index), 7))\n",
"ind_cent = np.zeros((len(common.index), 7))\n",
"\n",
"fractlist = ['3K', '5.4K', '12.2K', '24K', '78.4K', '110K', '195.5K']\n",
"\n",
"for protein, row in zip(common.index, range(0, len(common.index))):\n",
" for fraction, col in zip(fractlist, range(0, 7)):\n",
" tempa1 = ua1.at[protein, fraction]\n",
" tempb1 = ub1.at[protein, fraction]\n",
" tempc1 = uc1.at[protein, fraction]\n",
" tempa2 = ua1.at[protein, fraction]\n",
" tempb2 = ub1.at[protein, fraction]\n",
" tempc2 = uc1.at[protein, fraction]\n",
" un_cent[row, col] = np.mean([tempa1, tempb1, tempc1, tempa2, tempb2, tempc2])\n",
"\n",
"for protein, row in zip(common.index, range(0, len(common.index))):\n",
" for fraction, col in zip(fractlist, range(0, 7)):\n",
" tempa1 = ia2.at[protein, fraction]\n",
" tempb1 = ib2.at[protein, fraction]\n",
" tempc1 = ic2.at[protein, fraction]\n",
" tempa2 = ia2.at[protein, fraction]\n",
" tempb2 = ib2.at[protein, fraction]\n",
" tempc2 = ic2.at[protein, fraction]\n",
" ind_cent[row, col] = np.mean([tempa1, tempb1, tempc1, tempa2, tempb2, tempc2])"
]
},
{
"cell_type": "code",
"execution_count": 69,
"metadata": {},
"outputs": [],
"source": [
"# Calculate Euclidean distance from each replicate to its respective centroid\n",
"un_stats = np.zeros((len(common.index), 7))\n",
"ind_stats = np.zeros((len(common.index), 7))\n",
"\n",
"unlist = [ua1, ub1, uc1, ua2, ub2, uc2]\n",
"indlist = [ia1, ib1, ic1, ia2, ib2, ic2]\n",
"\n",
"for protein, row in zip(common.index, range(0, len(common.index))):\n",
" for df, col in zip(unlist, range(0,6)):\n",
" temp = df.loc[protein, '3K':'195.5K']\n",
" avg = un_cent[row]\n",
" un_stats[row, col] = euclidean(temp, avg)\n",
"\n",
"for protein, row in zip(common.index, range(0, len(common.index))):\n",
" for df, col in zip(indlist, range(0,6)):\n",
" temp = df.loc[protein, '3K':'195.5K']\n",
" avg = ind_cent[row]\n",
" ind_stats[row, col] = euclidean(temp, avg)"
]
},
{
"cell_type": "code",
"execution_count": 70,
"metadata": {},
"outputs": [],
"source": [
"# Calculate standard deviation of distances for each condition\n",
"un_stats[:, 6] = np.std(un_stats[:, 0:6], axis=1, ddof=1)\n",
"ind_stats[:, 6] = np.std(ind_stats[:, 0:6], axis=1, ddof=1)"
]
},
{
"cell_type": "code",
"execution_count": 71,
"metadata": {},
"outputs": [],
"source": [
"# Calculate t-statistic and p-value for each protein\n",
"# Using Benjamini-Hochberg correction for p-values\n",
"mvmt_stats = np.zeros((len(common.index), 5))\n",
"\n",
"# Column 0 will contain the Euclidean distance between the uninduced and induced centroids\n",
"for row in range(0, len(common.index)):\n",
" mvmt_stats[row, 0] = euclidean(un_cent[row], ind_cent[row])\n",
"\n",
"# Column 1 will contain the t-statistic\n",
"for row in range(0, len(common.index)):\n",
" unvar = (un_stats[row, 3])**2\n",
" indvar = (ind_stats[row, 3])**2\n",
" mvmt_stats[row, 1] = ((mvmt_stats[row, 0])/np.sqrt((unvar+indvar)/3))\n",
"\n",
"# Column 2 will contain the p-value\n",
"for row in range(0, len(common.index)):\n",
" tstat = mvmt_stats[row, 1]\n",
" unvar = (un_stats[row, 3])**2\n",
" indvar = (ind_stats[row, 3])**2\n",
" dof = (2*((unvar**2+(2*unvar*indvar)+indvar**2)/(unvar**2+indvar**2)))\n",
" mvmt_stats[row, 2] = t.sf(x=tstat, df=dof)\n",
"\n",
"# Column 3 will contain the rank value for each p-value\n",
"mvmt_stats[:, 3] = rankdata(mvmt_stats[:, 2])\n",
"\n",
"# Column 4 will contain the Benjamini-Hochberg critical values; using an FDR of 0.266 or 26.6%\n",
"mvmt_stats[:, 4] = (mvmt_stats[:, 3]/len(common.index))*0.266"
]
},
{
"cell_type": "code",
"execution_count": 72,
"metadata": {},
"outputs": [],
"source": [
"# Adding the p-value and critical value for each protein, then sorting by p-value; FDR of 26.6%\n",
"sigmvmt = common.loc[:, 'Gene_ua1':'Protein information_ua1']\n",
"sigmvmt['p-value'] = mvmt_stats[:, 2]\n",
"sigmvmt['CritValue'] = mvmt_stats[:, 4]\n",
"sigmvmt.sort_values(by='p-value', inplace=True)"
]
},
{
"cell_type": "code",
"execution_count": 73,
"metadata": {},
"outputs": [],
"source": [
"# Determining the maximum p-value < CritValue and subsetting out all p-values less than that cutoff; FDR of 26.6%\n",
"# Taking just the top 500\n",
"df = sigmvmt[sigmvmt['p-value'] < sigmvmt['CritValue']]\n",
"cutoff = df['p-value'].max()\n",
"moved = sigmvmt[sigmvmt['p-value'] <= cutoff]"
]
},
{
"cell_type": "code",
"execution_count": 74,
"metadata": {},
"outputs": [],
"source": [
"# Saving top 500 from combined movement analysis\n",
"moved500 = moved.iloc[0:500, :]\n",
"moved500.to_csv('20190528_CombinedMvmtAnalysis_top500_26.6FDR.csv', sep=',')\n",
"moved500.reset_index(inplace=True)\n",
"df = moved500.iloc[:, 0:1]\n",
"df.to_csv('20190528_CombinedMvmtAnalysis_top500_IDs.csv', sep=',', header=False, index=False)"
]
},
{
"cell_type": "code",
"execution_count": 75,
"metadata": {},
"outputs": [],
"source": [
"# Comparing the 1.5% FDR list of WT2 and the 2% FDR list of WT1, only looking at top 1,000 proteins from each\n",
"wt1_top1000 = wt1_moved.iloc[0:1000, :]\n",
"wt2_top1000 = wt2_moved.iloc[0:1000, :]\n",
"moved_common = pd.concat([wt1_top1000, wt2_top1000], join='inner', axis=1)\n",
"moved_common.to_csv('20190605_ComparativeMvmtAnalysis.csv', sep=',')"
]
},
{
"cell_type": "code",
"execution_count": 76,
"metadata": {},
"outputs": [],
"source": [
"# Making IDs only files for moved_common\n",
"moved_common.reset_index(inplace=True)\n",
"df = moved_common.iloc[:, 0:1]\n",
"df.to_csv('20190605_ComparativeMvmtAnalysis_IDs.csv', sep=',', header=False, index=False)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## STRING Analyses for Centroid Distance Analysis\n",
"\n",
"[Top 500 proteins from combined distance analysis](https://version-11-0.string-db.org/cgi/network.pl?networkId=YVBEOmahI0G8)<br>\n",
"[Common moved proteins between the biological replicates](https://version-11-0.string-db.org/cgi/network.pl?networkId=jz0AgR89gEGY)"
]
}
],
"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.7.4"
}
},
"nbformat": 4,
"nbformat_minor": 4
}