{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# What is about ? \n\n\nHere we analyse rna which are coexpressed with CD32 gene, both for its rna and the protein\n\nFor the begining we will mainly use the Kaggle-NIPS 2022 competition dataset.\n\n\n\n---------------------------------------------------\n\nThe project:\n\nFar goals:\n\n    1) create \"full\" list of rna related to the selected protein and its correposponding RNA\n    2) compare them and try to relate the   diffence to post-transcriptional, post-translational regulation\n\nGoals:\n\n    1) calculate feature importances by several methods and compare them\n    1.1) univariable importances - correlations, chi2, mutual info, etc..\n    1.2) model feature importances - Ridge,Lasso, Permutation FI, Boostings, Shap ...\n    1.3) for each importance - analyse stability with various respects - subsample, model paramers\n    1.4) try to understand - how much these importances add to basic ones - like correlations\n    1.5) do ML importances give something new comparing to univariable ? May be yes - if one have \"synergy\" , also due to specificity of data - many \"dropouts\" - \"collective\" methods like - ML models may work better. \n    1.6) ML can be faster than mutual info e.g. - analyse computational resources\n    1.7) find \"all\" features meaningfully related to tartges\n    1.8) how reliable are ML methods - they can choose features which actually are NOT related to target - how to overcome it ? \n    1.9) how to select the actual number of important methods features in ML methods ? how many to take ? \n\n    2) Try to find interpretation of obtained importances\n    2.1) cluster them - may be clusters - will give clear picture\n    2.2) enrichment analysis\n    2.3) look on umap for localization ","metadata":{"editable":false}},{"cell_type":"markdown","source":"# Preparations","metadata":{"editable":false}},{"cell_type":"code","source":"# This Python 3 environment comes with many helpful analytics libraries installed\n# It is defined by the kaggle/python Docker image: https://github.com/kaggle/docker-python\n# For example, here's several helpful packages to load\n\nimport matplotlib.pyplot as plt\nimport seaborn as sns\nimport time \nt0start = time.time()\nimport gc\n\nimport numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\n\n# Input data files are available in the read-only \"../input/\" directory\n# For example, running this (by clicking run or pressing Shift+Enter) will list all files under the input directory\n\nimport os\nfor dirname, _, filenames in os.walk('/kaggle/input'):\n    for filename in filenames:\n        if 'research-project-01' not in os.path.join(dirname, filename): \n            print(os.path.join(dirname, filename))\n\n# You can write up to 20GB to the current directory (/kaggle/working/) that gets preserved as output when you create a version using \"Save & Run All\" \n# You can also write temporary files to /kaggle/temp/, but they won't be saved outside of the current session","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","editable":false,"execution":{"iopub.status.busy":"2023-04-18T12:55:34.533060Z","iopub.execute_input":"2023-04-18T12:55:34.534237Z","iopub.status.idle":"2023-04-18T12:55:35.734189Z","shell.execute_reply.started":"2023-04-18T12:55:34.534184Z","shell.execute_reply":"2023-04-18T12:55:35.732736Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Load data","metadata":{"editable":false}},{"cell_type":"code","source":"str_data_inf = 'NIPS22'","metadata":{"editable":false,"execution":{"iopub.status.busy":"2023-04-18T12:55:35.736635Z","iopub.execute_input":"2023-04-18T12:55:35.737204Z","iopub.status.idle":"2023-04-18T12:55:35.743095Z","shell.execute_reply.started":"2023-04-18T12:55:35.737155Z","shell.execute_reply":"2023-04-18T12:55:35.741733Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\n\n# Here we use more general setup which allows to with many datasets - from the Notebook: https://www.kaggle.com/code/alexandervc/many-citeseq-02-just-load-full-data\n# Presently we consider only one dataset so it is overkill, but nevertheless \n\ndict_df_prot = {}\ndict_df_rna = {}\ndict_df_meta = {}\n\nN_rows2take = int( 1e5 ) # Work only with first rows of the data - to speed up, avoid RAM crashes , etc... \nwork_with_rna_data = 'Yes'# If No - we will not store/work with protein data - only rna -  For some fast analysis we need only protein data, which much more smaller \nwork_with_protein_data = 'Yes' \ndict_datasets2consider = {}\n# dict_datasets2consider['Kaggle2302'] = 'Yes'\ndict_datasets2consider['NIPS2022'] = 'Yes'\n\n#  %%time\nimport gc\nif dict_datasets2consider['NIPS2022'] == 'Yes':\n\n    key4dict = 'NIPS2022'\n    \n    DATA_DIR = \"/kaggle/input/open-problems-multimodal/\"\n    if work_with_protein_data == 'Yes':\n        FP_CITE_TRAIN_TARGETS = os.path.join(DATA_DIR,\"train_cite_targets.h5\")\n        df = pd.read_hdf(FP_CITE_TRAIN_TARGETS, start=0, stop=N_rows2take, )\n        print( [t for t in df.columns if 'CD45' in t.upper() ], [t for t in df.columns if 'CD53' in t.upper() ],  )\n        dict_df_prot[key4dict] = df\n        print('Protein data shape:', df.shape)\n        display(df.head(2)) \n\n    # %%time\n    if work_with_rna_data == 'Yes':\n        filename_rna_data = '/kaggle/input/open-problems-multimodal/train_cite_inputs.h5'\n        df = pd.read_hdf(filename_rna_data, start=0, stop=N_rows2take,)\n        l = [t.split('_')[1] for t in df.columns]; df.columns = l \n        index_save = df.index\n        print('RNA data shape:', df.shape)\n        dict_df_rna[key4dict] = df\n        display(df.head(2)) \n\n\n    fn = '/kaggle/input/open-problems-multimodal/metadata.csv'\n    df = pd.read_csv(fn, index_col = 0 );\n    mask = df.index.isin(index_save); df = df[mask]\n    print('Meta data shape:', df.shape)\n    dict_df_meta[key4dict] = df\n    display(df.head(2)) \n\n    gc.collect() \n\n\nfor k in dict_df_meta:\n    df_prot = dict_df_prot[k]\n    df_rna = dict_df_rna[k]\n    df_meta = dict_df_meta[k]\n    print(df_prot.shape, df_rna.shape, df_meta.shape )","metadata":{"editable":false,"execution":{"iopub.status.busy":"2023-04-18T12:55:35.745176Z","iopub.execute_input":"2023-04-18T12:55:35.746064Z","iopub.status.idle":"2023-04-18T12:56:31.597199Z","shell.execute_reply.started":"2023-04-18T12:55:35.746006Z","shell.execute_reply":"2023-04-18T12:56:31.595870Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"protein_name =  'CD32'# 'CD161'#'CD8a' #  'CD45RA'# 'CD45RO'#  'HLA.DR'\nrna_name = 'FCGR2A' #  dict_CD_proteinName_to_CD_rnaName[protein_name]\nprint(protein_name, rna_name )","metadata":{"editable":false,"execution":{"iopub.status.busy":"2023-04-18T12:56:31.599793Z","iopub.execute_input":"2023-04-18T12:56:31.600235Z","iopub.status.idle":"2023-04-18T12:56:31.606769Z","shell.execute_reply.started":"2023-04-18T12:56:31.600194Z","shell.execute_reply":"2023-04-18T12:56:31.604997Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# First look","metadata":{"editable":false}},{"cell_type":"code","source":"print('Protein:', protein_name, ' rna: ', rna_name )\nx = df_prot[protein_name]\ny = df_rna[rna_name]\nd = pd.DataFrame(); d['Prot'] = x ; d['Rna'] = y\nprint('Count zeros:',  dict( (d==0).sum(axis=0)  ), 'Percents of zeros: ', dict( np.round( (d==0).sum(axis=0)/len(d)*100, 1)  ),  )\nprint('Correlations:')\ndisplay(d.corr() )\nmask = (d.iloc[:,0] != 0)& (d.iloc[:,1] != 0 )\nprint('Correlations on part without zeros (\"dropout\"):')\ndisplay(d[mask].corr() )\n\ndisplay(d.describe() )\nfig = plt.figure ( figsize = ( 20,6 ) )\nplt.suptitle(str_data_inf , fontsize = 20 )\nfig.add_subplot(1,2,1)\nplt.hist( d['Prot'], bins = 100 ) \nplt.title('Protein '+protein_name, fontsize = 20 )\nfig.add_subplot(1,2,2)\nplt.hist( d['Rna'] , bins = 100 ) \nplt.title('RNA '+rna_name, fontsize = 20 )\nplt.show()\ns = (d['Prot'] > 0.5).sum(axis=0)\nprint('Protein > 0.5:',  ( (d['Prot'] > 0.5).sum(axis=0)  ), 'Percent Protein > 0.5: ', ( np.round( s /len(d)*100, 1)  ),  )\ns = (d['Rna'] > 0).sum(axis=0)\nprint('Rna > 0 :',  s, 'Percent Rna > 0: ', ( np.round( s /len(d)*100, 1)  ),  )\n\n\nfig = plt.figure ( figsize = ( 10,10 ) )\nsns.scatterplot(x = d['Prot'] , y = d['Rna'] )\nplt.title(str_data_inf + ' ' + protein_name  + ' ' +  rna_name , fontsize = 20   )\nplt.show()\n\nfig = plt.figure ( figsize = ( 10,10 ) )\n# sns.scatterplot(x = d['Prot'] , y = d['Rna'] )\nplt.hist2d(x = d['Prot'] , y = d['Rna'], bins = 20 )\nplt.title(str_data_inf + ' ' + protein_name  + ' ' +  rna_name , fontsize = 20   )\nplt.show()\n\nfig = plt.figure ( figsize = ( 10,10 ) )\n# sns.scatterplot(x = d['Prot'] , y = d['Rna'] )\nmask =  d['Rna'] > 0 \nplt.hist2d(x = d['Prot'][mask] , y = d['Rna'][mask] , bins = 20)\n\nplt.title(str_data_inf + ' ' + protein_name  + ' ' +  rna_name , fontsize = 20   )\nplt.show()\n","metadata":{"editable":false,"execution":{"iopub.status.busy":"2023-04-18T12:56:31.608977Z","iopub.execute_input":"2023-04-18T12:56:31.609467Z","iopub.status.idle":"2023-04-18T12:56:33.153394Z","shell.execute_reply.started":"2023-04-18T12:56:31.609422Z","shell.execute_reply":"2023-04-18T12:56:33.152262Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Function to analyse importances","metadata":{"editable":false}},{"cell_type":"code","source":"# Based on: \n# https://www.kaggle.com/code/alexandervc/mmscel-feature-importance-03-univariable?scriptVersionId=119140380&cellId=52\n    \ndef importances_stat(df_importances, duplicate_plots_restricted_to_N = 500, \n                     list_topN_to_take_for_df_stat = [100,150, 200, 300, 1000],\n                     make_plots = 'All'\n                    ):\n\n    if make_plots == 'All':\n        n_x_subplots = 2; c = 0 \n        for col in df_importances.columns:\n            if c % n_x_subplots == 0:\n                if c > 0:\n                    plt.show()\n                fig = plt.figure(figsize = (20,5) ); c = 0    \n            c += 1; fig.add_subplot(1,n_x_subplots ,c)\n            if 'Corr' in col:\n                m = df_importances[col] < 0.999 # Exclude correlation equal to 1 since typically it is correlation with itself\n            else:\n                m = pd.Series(data = True, index = df_importances.index )\n\n            plt.hist(df_importances[col][m], bins = 100)\n            plt.title(col)\n        plt.show()\n\n    display(df_importances.describe() )    \n    col = df_importances.columns[0]\n    df_importances.sort_values(col, ascending = False, key = abs ).head(20)\n\n\n    df_stat = pd.DataFrame()\n    N = df_importances.shape[1]\n    for i in range(N):\n        for j in range(N):\n            if j<=i: continue\n            col1 = df_importances.columns[i];         col2 = df_importances.columns[j];\n            l = []\n            if 'LGB' not in col1:\n                v1 = df_importances[col1].sort_values(ascending = False, key = abs )\n            else:\n                v1 = df_importances[col1].sort_values(ascending = False)\n            if 'LGB' not in col2:\n                v2 = df_importances[col2].sort_values(ascending = False, key = abs )\n            else:\n                v2 = df_importances[col2].sort_values(ascending = False)\n\n            for k in list_topN_to_take_for_df_stat:# range(len(v1)): # [10,25,50, 75, 100, 1000, 2000, 4000]:\n                s = set(v1.index[:k] ) & set(v2.index[:k] ) \n                l.append(len(s) )\n                # print(k, len(s))   \n\n            print(col1,col2, len(l))\n            for  ii,k in enumerate(list_topN_to_take_for_df_stat):\n                df_stat.loc[col1+' VS '+col2, 'Top'+str(k)+ ' Intersection'] = l[ii]\n\n    display(df_stat )\n\n    N = df_importances.shape[1]\n    for i in range(N):\n        for j in range(N):\n            if j<=i: continue\n            col1 = df_importances.columns[i];         col2 = df_importances.columns[j];\n\n            if 'LGB' not in col1:\n                df_importances_ranks = (-df_importances.sort_values(col1, ascending = False, key = abs).abs()).rank()\n            else:\n                df_importances_ranks = (-df_importances.sort_values(col1, ascending = False)).rank()\n\n            v1 = df_importances_ranks[col1] - df_importances_ranks[col2]\n\n            for  ii,k in enumerate(list_topN_to_take_for_df_stat):\n                df_stat.loc[col1+' VS '+col2, 'Max Rank Diffr Top'+str(k)] = np.max( np.abs( v1.values[:k] ) )\n            \n\n    display(df_stat)\n\n    if make_plots == 'All':\n        #%%time\n        N = df_importances.shape[1]\n        for i in range(N):\n            for j in range(N):\n                if j<=i: continue\n                col1 = df_importances.columns[i];         col2 = df_importances.columns[j];\n                print(col1,col2)\n                l = []\n                if 'LGB' not in col1:\n                    v1 = df_importances[col1].sort_values(ascending = False, key = abs )\n                else:\n                    v1 = df_importances[col1].sort_values(ascending = False)\n                if 'LGB' not in col2:\n                    v2 = df_importances[col2].sort_values(ascending = False, key = abs )\n                else:\n                    v2 = df_importances[col2].sort_values(ascending = False)\n\n                for k in range(len(v1)): # [10,25,50, 75, 100, 1000, 2000, 4000]:\n                    s = set(v1.index[:k] ) & set(v2.index[:k] ) \n                    l.append(len(s) )\n                    # print(k, len(s))       \n\n                fig = plt.figure(figsize = (20,6))\n                plt.title(col1+' ' + col2, fontsize = 20 )\n                plt.plot(l ,label = 'original data - counts intersection' )\n\n                y = np.array(l)\n                x = np.arange(len(y))\n                p = np.polyfit(x,y,2)\n                print(p)\n                plt.plot(np.polyval(p,x),label = 'parabolic approximation' )\n\n                y = np.array(l)\n                x = np.arange(len(y))\n                p = np.polyfit(x[:50],y[:50],1)\n                print(p)\n                plt.plot(np.polyval(p,x),label = 'linear approximation  based on 0-50 points' )\n                \n                \n                i0 =  int( len(x)/2 )\n                l_tmp = np.arange(0,i0, 100)\n                y = np.array(l)\n                x = np.arange(len(y))\n                #print(l_tmp,  x[l_tmp], y[l_tmp] )\n                p = np.polyfit(x[l_tmp],y[l_tmp],1)\n                print('linear for half size:',p)\n                plt.plot(np.polyval(p,x),label = 'linear approximation  based on half of the cells' )\n\n                plt.legend(fontsize = 20 )\n                plt.grid()\n                plt.show()            \n\n\n                if (duplicate_plots_restricted_to_N is not None ) and (duplicate_plots_restricted_to_N > 0):\n                    N_loc = duplicate_plots_restricted_to_N\n                    fig = plt.figure(figsize = (20,10) )\n                    plt.title(col1+' ' + col2 , fontsize = 20)\n                    plt.plot(l ,label = 'original data - counts intersection' )\n\n                    y = np.array(l)\n                    x = np.arange(len(y))\n                    p = np.polyfit(x,y,2)\n                    print(p)\n                    plt.plot(np.polyval(p,x),label = 'parabolic approximation' )\n\n                    y = np.array(l)\n                    x = np.arange(len(y))\n                    p = np.polyfit(x[:50],y[:50],1)\n                    print(p)\n                    plt.plot(np.polyval(p,x),label = 'linear approximation based on 0-50 points' )\n\n                    y = np.array(l)\n                    x = np.arange(len(y))\n                    p = np.polyfit(x[:150],y[:150],1)\n                    print(p)\n                    plt.plot(np.polyval(p,x),label = 'linear approximation  based on 0-150 points' )\n\n                    plt.xlim([0,N_loc])\n                    plt.ylim([0,N_loc/2])\n\n                    plt.legend(fontsize = 20 )\n                    plt.grid()\n                    plt.show()        \n\n                ##################################################################33\n                l2 = []\n                for k in range(10,len(l)):\n                    l2.append(l[k]/k)\n\n                fig = plt.figure(figsize = (20,10) )\n                plt.title(col1+' ' + col2 , fontsize = 20)\n                plt.plot(l2 ,label = 'intersection topK / K' )\n\n                plt.legend(fontsize = 20 )\n                plt.grid()\n                plt.show()                \n\n                if (duplicate_plots_restricted_to_N is not None ) and (duplicate_plots_restricted_to_N > 0):\n                    N_loc = duplicate_plots_restricted_to_N\n                    fig = plt.figure(figsize = (20,10) )\n                    plt.title(col1+' ' + col2 , fontsize = 20)\n                    plt.plot(l2[:N_loc] ,label = 'intersection topK / K' )\n\n                    plt.xlim([0,N_loc])\n\n                    plt.legend(fontsize = 20 )\n                    plt.grid()\n                    plt.show()                        \n\n        N = df_importances.shape[1]\n        for i in range(N):\n            for j in range(N):\n                if j<=i: continue\n                col1 = df_importances.columns[i];         col2 = df_importances.columns[j];\n\n                if 'LGB' not in col1:\n                    df_importances_ranks = (-df_importances.sort_values(col1, ascending = False, key = abs).abs()).rank()\n                else:\n                    df_importances_ranks = (-df_importances.sort_values(col1, ascending = False)).rank()\n\n                v1 = df_importances_ranks[col1] - df_importances_ranks[col2]\n\n                fig = plt.figure(figsize = (20,6))\n                plt.title('Ranks differences: ' + col1+' ' + col2, fontsize = 20 )\n                plt.plot(v1.values ) \n                plt.grid()\n                plt.show()\n\n                if (duplicate_plots_restricted_to_N is not None ) and (duplicate_plots_restricted_to_N > 0):\n                    N_loc = duplicate_plots_restricted_to_N\n                    fig = plt.figure(figsize = (20,6))\n                    plt.title('Ranks differences: ' + col1+' ' + col2, fontsize = 20 )\n                    plt.plot(v1.values[:N_loc] ) \n                    plt.grid()\n                    plt.show()                \n\n    df_stat.to_csv('stat_importances_'+str_data_inf+str(protein_name).replace('.','-')+'_'+ str(rna_name).replace('.','-') +'.csv')            \n    return df_stat","metadata":{"editable":false,"execution":{"iopub.status.busy":"2023-04-18T12:56:33.155208Z","iopub.execute_input":"2023-04-18T12:56:33.155943Z","iopub.status.idle":"2023-04-18T12:56:33.203924Z","shell.execute_reply.started":"2023-04-18T12:56:33.155892Z","shell.execute_reply":"2023-04-18T12:56:33.202493Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Calculate importances by Pearson correlation ","metadata":{"editable":false}},{"cell_type":"code","source":"%%time\nfrom sklearn.feature_selection import r_regression,  f_regression, mutual_info_regression\n\ndf_importances = pd.DataFrame(index = df_rna.columns) \n\nvec_corr = r_regression(df_rna, x )\ndf_importances[protein_name +' Prot Corr'] = vec_corr\n\nvec_corr = r_regression(df_rna, y )\ndf_importances[rna_name +' RNA Corr'] = vec_corr","metadata":{"editable":false,"execution":{"iopub.status.busy":"2023-04-18T12:56:33.205987Z","iopub.execute_input":"2023-04-18T12:56:33.206481Z","iopub.status.idle":"2023-04-18T12:56:49.681222Z","shell.execute_reply.started":"2023-04-18T12:56:33.206430Z","shell.execute_reply":"2023-04-18T12:56:49.676882Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print('Features counts', len(df_importances ) )\nprint('Percent of features abs-correlated higher than 0.1:')\nprint( (df_importances.abs() > 0.1 ).sum(axis = 0)/ len(df_importances) * 100 )\nprint()\nprint('Count of features abs-correlated higher than 0.1:')\nprint( (df_importances.abs() > 0.1 ).sum(axis = 0))# / len(df_importances) * 100 )\n\ncol = df_importances.columns[0]\nd1 = df_importances.sort_values(col,  ascending = False, key = abs ).reset_index()\ncol = df_importances.columns[1]\nd2 = df_importances.sort_values(col,  ascending = False, key = abs ).reset_index()\ndisplay( pd.concat([d1,d2], axis = 1).head(20) )\npd.concat([d1,d2], axis = 1).tail(3)","metadata":{"editable":false,"execution":{"iopub.status.busy":"2023-04-18T12:56:49.688304Z","iopub.execute_input":"2023-04-18T12:56:49.689985Z","iopub.status.idle":"2023-04-18T12:56:49.774525Z","shell.execute_reply.started":"2023-04-18T12:56:49.689884Z","shell.execute_reply":"2023-04-18T12:56:49.773018Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_stat = importances_stat( df_importances, duplicate_plots_restricted_to_N = None )\ndisplay(df_stat)","metadata":{"editable":false,"execution":{"iopub.status.busy":"2023-04-18T12:56:49.776431Z","iopub.execute_input":"2023-04-18T12:56:49.776938Z","iopub.status.idle":"2023-04-18T12:58:26.201773Z","shell.execute_reply.started":"2023-04-18T12:56:49.776903Z","shell.execute_reply":"2023-04-18T12:58:26.200445Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Pearson correlation excluding cells with zero values of protein or rna (dropout)","metadata":{"editable":false}},{"cell_type":"code","source":"%%time\nfrom sklearn.feature_selection import r_regression,  f_regression, mutual_info_regression\n\nmask = (x!=0) & (y!=0)\n\ndf_importances_loc = pd.DataFrame(index = df_rna.columns) \n\nvec_corr = r_regression(df_rna[mask], x[mask] )\ndf_importances[protein_name +' Prot Corr NoZeros (Dropout)'] = vec_corr\ndf_importances_loc[protein_name +' Prot Corr NoZeros (Dropout)'] = vec_corr\n\n\nvec_corr = r_regression(df_rna[mask], y[mask] )\ndf_importances[rna_name +' RNA Corr  NoZeros (Dropout)'] = vec_corr\ndf_importances_loc[rna_name +' RNA Corr  NoZeros (Dropout)'] = vec_corr\n\n\nprint('Features counts', len(df_importances ) )\nprint('Percent of features abs-correlated higher than 0.1:')\nprint( (df_importances.abs() > 0.1 ).sum(axis = 0)/ len(df_importances) * 100 )\nprint()\nprint('Count of features abs-correlated higher than 0.1:')\nprint( (df_importances.abs() > 0.1 ).sum(axis = 0))# / len(df_importances) * 100 )","metadata":{"editable":false,"execution":{"iopub.status.busy":"2023-04-18T12:58:26.207345Z","iopub.execute_input":"2023-04-18T12:58:26.207971Z","iopub.status.idle":"2023-04-18T12:58:31.646250Z","shell.execute_reply.started":"2023-04-18T12:58:26.207926Z","shell.execute_reply":"2023-04-18T12:58:31.645075Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_stat = importances_stat( df_importances_loc , duplicate_plots_restricted_to_N = None)\ndisplay(df_stat)","metadata":{"editable":false,"execution":{"iopub.status.busy":"2023-04-18T12:58:31.647890Z","iopub.execute_input":"2023-04-18T12:58:31.648339Z","iopub.status.idle":"2023-04-18T13:00:08.805172Z","shell.execute_reply.started":"2023-04-18T12:58:31.648293Z","shell.execute_reply":"2023-04-18T13:00:08.803927Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Calculate importances by Spearman correlation ","metadata":{"editable":false}},{"cell_type":"code","source":"%%time\n\ndf_importances_spearman = pd.DataFrame(index = df_rna.columns) \n\nvec_corr = df_rna.corrwith(x, method=\"spearman\")\ndf_importances_spearman[protein_name +' Prot Spearman Corr'] = vec_corr\n\nvec_corr = df_rna.corrwith(y, method=\"spearman\")\ndf_importances_spearman[rna_name +' RNA Spearman Corr'] = vec_corr\n\ngc.collect()","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:00:08.808044Z","iopub.execute_input":"2023-04-18T13:00:08.809096Z","iopub.status.idle":"2023-04-18T13:12:36.747059Z","shell.execute_reply.started":"2023-04-18T13:00:08.809042Z","shell.execute_reply":"2023-04-18T13:12:36.745628Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print('Features counts', len(df_importances_spearman ) )\nprint('Percent of features abs-correlated higher than 0.1:')\nprint( (df_importances_spearman.abs() > 0.1 ).sum(axis = 0)/ len(df_importances_spearman) * 100 )\nprint()\nprint('Count of features abs-correlated higher than 0.1:')\nprint( (df_importances_spearman.abs() > 0.1 ).sum(axis = 0))# / len(df_importances) * 100 )\n\ncol = df_importances_spearman.columns[0]\nd1 = df_importances_spearman.sort_values(col,  ascending = False, key = abs ).reset_index()\ncol = df_importances_spearman.columns[1]\nd2 = df_importances_spearman.sort_values(col,  ascending = False, key = abs ).reset_index()\ndisplay( pd.concat([d1,d2], axis = 1).head(20) )\npd.concat([d1,d2], axis = 1).tail(3)","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:12:36.748878Z","iopub.execute_input":"2023-04-18T13:12:36.750212Z","iopub.status.idle":"2023-04-18T13:12:36.806271Z","shell.execute_reply.started":"2023-04-18T13:12:36.750159Z","shell.execute_reply":"2023-04-18T13:12:36.804945Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_stat_spearman = importances_stat( df_importances_spearman, duplicate_plots_restricted_to_N = None )\ndisplay(df_importances_spearman)","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:12:36.807989Z","iopub.execute_input":"2023-04-18T13:12:36.808360Z","iopub.status.idle":"2023-04-18T13:14:18.432063Z","shell.execute_reply.started":"2023-04-18T13:12:36.808322Z","shell.execute_reply":"2023-04-18T13:14:18.430736Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Spearman correlation excluding cells with zero values of protein or rna (dropout)","metadata":{}},{"cell_type":"code","source":"%%time\n\nmask = (x!=0) & (y!=0)\n\nvec_corr = df_rna[mask].corrwith(x, method=\"spearman\")\ndf_importances_spearman[protein_name +' Prot Spearman Corr NoZeros (Dropout)'] = vec_corr\n\n\nvec_corr = df_rna.corrwith(y, method=\"spearman\")\ndf_importances_spearman[rna_name +' RNA Spearman Corr  NoZeros (Dropout)'] = vec_corr\n\n\nprint('Features counts', len(df_importances ) )\nprint('Percent of features abs-correlated higher than 0.1:')\nprint( (df_importances_spearman.abs() > 0.1 ).sum(axis = 0)/ len(df_importances) * 100 )\nprint()\nprint('Count of features abs-correlated higher than 0.1:')\nprint( (df_importances_spearman.abs() > 0.1 ).sum(axis = 0))# / len(df_importances) * 100 )","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:14:18.433886Z","iopub.execute_input":"2023-04-18T13:14:18.434232Z","iopub.status.idle":"2023-04-18T13:23:36.907196Z","shell.execute_reply.started":"2023-04-18T13:14:18.434198Z","shell.execute_reply":"2023-04-18T13:23:36.905755Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_stat_spearman = importances_stat( df_importances_spearman, duplicate_plots_restricted_to_N = None)\ndisplay(df_stat_spearman)","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:23:36.908783Z","iopub.execute_input":"2023-04-18T13:23:36.909504Z","iopub.status.idle":"2023-04-18T13:34:11.897097Z","shell.execute_reply.started":"2023-04-18T13:23:36.909465Z","shell.execute_reply":"2023-04-18T13:34:11.895837Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_stat = importances_stat( df_importances_loc , duplicate_plots_restricted_to_N = None)\ndisplay(df_stat)","metadata":{"editable":false,"execution":{"iopub.status.busy":"2023-04-18T13:34:11.898875Z","iopub.execute_input":"2023-04-18T13:34:11.899676Z","iopub.status.idle":"2023-04-18T13:35:47.776245Z","shell.execute_reply.started":"2023-04-18T13:34:11.899603Z","shell.execute_reply":"2023-04-18T13:35:47.774998Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Show stat on all importances ","metadata":{"editable":false}},{"cell_type":"code","source":"df_stat = importances_stat( df_importances , duplicate_plots_restricted_to_N = 400, make_plots = None)\n#display(df_stat)","metadata":{"editable":false,"execution":{"iopub.status.busy":"2023-04-18T13:35:47.778139Z","iopub.execute_input":"2023-04-18T13:35:47.778889Z","iopub.status.idle":"2023-04-18T13:35:47.997589Z","shell.execute_reply.started":"2023-04-18T13:35:47.778830Z","shell.execute_reply":"2023-04-18T13:35:47.996479Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print('%.1f seconds passed total '%(time.time()-t0start) )","metadata":{"editable":false,"execution":{"iopub.status.busy":"2023-04-18T13:35:47.999099Z","iopub.execute_input":"2023-04-18T13:35:47.999554Z","iopub.status.idle":"2023-04-18T13:35:48.006247Z","shell.execute_reply.started":"2023-04-18T13:35:47.999520Z","shell.execute_reply":"2023-04-18T13:35:48.004923Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Check non-overlapping genes","metadata":{}},{"cell_type":"code","source":"!pip install gseapy","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:35:48.008091Z","iopub.execute_input":"2023-04-18T13:35:48.008598Z","iopub.status.idle":"2023-04-18T13:36:02.055598Z","shell.execute_reply.started":"2023-04-18T13:35:48.008538Z","shell.execute_reply":"2023-04-18T13:36:02.054531Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import gseapy as gp\n\ndef get_ntop(df, col, n=100):\n    return df[col].abs().sort_values(ascending=False)[:n].index","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:36:02.057944Z","iopub.execute_input":"2023-04-18T13:36:02.058832Z","iopub.status.idle":"2023-04-18T13:36:02.086753Z","shell.execute_reply.started":"2023-04-18T13:36:02.058777Z","shell.execute_reply":"2023-04-18T13:36:02.085462Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"top100_rna, top100_prot = get_ntop(df_importances_spearman, \"FCGR2A RNA Spearman Corr\"), \\\n                            get_ntop(df_importances_spearman, \"CD32 Prot Spearman Corr\")\ntop100_rna_non_overlap, top100_prot_non_overlap = list(set(top100_rna).difference(set(top100_prot))), \\\n                                                    list(set(top100_prot).difference(set(top100_rna)))\nprint(\"Non-overlapping RNA: {}\".format(len(top100_rna_non_overlap)))\nprint(\"Non-overlapping proteins: {}\".format(len(top100_prot_non_overlap)))","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:36:02.088476Z","iopub.execute_input":"2023-04-18T13:36:02.089288Z","iopub.status.idle":"2023-04-18T13:36:02.113282Z","shell.execute_reply.started":"2023-04-18T13:36:02.089236Z","shell.execute_reply":"2023-04-18T13:36:02.112336Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"gene_sets = [\"GO_Biological_Process_2021\", \"Reactome_2022\", \"WikiPathway_2021_Human\", \"KEGG_2021_Human\"]","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:36:02.114794Z","iopub.execute_input":"2023-04-18T13:36:02.115703Z","iopub.status.idle":"2023-04-18T13:36:02.120261Z","shell.execute_reply.started":"2023-04-18T13:36:02.115638Z","shell.execute_reply":"2023-04-18T13:36:02.119216Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## RNA Enrichment","metadata":{}},{"cell_type":"markdown","source":"Overcome unreachable Enrichr","metadata":{}},{"cell_type":"code","source":"!curl http://www.gsea-msigdb.org/gsea/msigdb/download_file.jsp?filePath=/msigdb/release/2023.1.Hs/msigdb_v2023.1.Hs_files_to_download_locally.zip","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:36:02.121879Z","iopub.execute_input":"2023-04-18T13:36:02.122555Z","iopub.status.idle":"2023-04-18T13:36:03.525871Z","shell.execute_reply.started":"2023-04-18T13:36:02.122519Z","shell.execute_reply":"2023-04-18T13:36:03.524340Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"os.listdir(\"/kaggle/working\")","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:36:03.528451Z","iopub.execute_input":"2023-04-18T13:36:03.529003Z","iopub.status.idle":"2023-04-18T13:36:03.539377Z","shell.execute_reply.started":"2023-04-18T13:36:03.528944Z","shell.execute_reply":"2023-04-18T13:36:03.538290Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"try:\n    rna_enrich = gp.enrichr(top100_rna_non_overlap, gene_sets=gene_sets, organism=\"Human\", cutoff=0.01)\nexcept:\n    rna_enrich = gp.enrichr(top100_rna_non_overlap, gene_sets=gene_sets, organism=\"Human\", cutoff=0.01)\ndisplay(rna_enrich.results.sort_values(\"Adjusted P-value\").head(10))","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:40:41.520432Z","iopub.execute_input":"2023-04-18T13:40:41.520863Z","iopub.status.idle":"2023-04-18T13:40:45.201284Z","shell.execute_reply.started":"2023-04-18T13:40:41.520824Z","shell.execute_reply":"2023-04-18T13:40:45.199920Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Protein Enrichment","metadata":{}},{"cell_type":"code","source":"try:\n    prot_enrich = gp.enrichr(top100_prot_non_overlap, gene_sets=gene_sets, organism=\"Human\", cutoff=0.01)\nexcept:\n    prot_enrich = gp.enrichr(top100_prot_non_overlap, gene_sets=gene_sets, organism=\"Human\", cutoff=0.01)\ndisplay(prot_enrich.results.sort_values(\"Adjusted P-value\").head(10))","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:40:49.672072Z","iopub.execute_input":"2023-04-18T13:40:49.672586Z","iopub.status.idle":"2023-04-18T13:40:54.980164Z","shell.execute_reply.started":"2023-04-18T13:40:49.672534Z","shell.execute_reply":"2023-04-18T13:40:54.978975Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Clustering of non-overlap RNAs","metadata":{}},{"cell_type":"code","source":"from sklearn.cluster import AgglomerativeClustering","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:40:58.826327Z","iopub.execute_input":"2023-04-18T13:40:58.827292Z","iopub.status.idle":"2023-04-18T13:40:59.005406Z","shell.execute_reply.started":"2023-04-18T13:40:58.827250Z","shell.execute_reply":"2023-04-18T13:40:59.004352Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def get_cluster_table(labels, clustering):\n    return pd.DataFrame(clustering.labels_, index=labels, columns=[\"cluster\"]).sort_values(\"cluster\")","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:41:00.027711Z","iopub.execute_input":"2023-04-18T13:41:00.028557Z","iopub.status.idle":"2023-04-18T13:41:00.035760Z","shell.execute_reply.started":"2023-04-18T13:41:00.028503Z","shell.execute_reply":"2023-04-18T13:41:00.034663Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sns.clustermap(df_rna[top100_rna_non_overlap].corr()) ","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:41:01.377253Z","iopub.execute_input":"2023-04-18T13:41:01.378408Z","iopub.status.idle":"2023-04-18T13:41:02.279037Z","shell.execute_reply.started":"2023-04-18T13:41:01.378352Z","shell.execute_reply":"2023-04-18T13:41:02.277696Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"clustering_rna = AgglomerativeClustering(linkage=\"average\", n_clusters=2, \n                                     compute_full_tree=True)\nclustering_rna.fit_predict(df_rna[top100_rna_non_overlap].corr())\ntop100_rna_non_overlap_clusters = get_cluster_table(df_rna[top100_rna_non_overlap].corr().index, clustering_rna)","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:41:03.785109Z","iopub.execute_input":"2023-04-18T13:41:03.785513Z","iopub.status.idle":"2023-04-18T13:41:03.978802Z","shell.execute_reply.started":"2023-04-18T13:41:03.785476Z","shell.execute_reply":"2023-04-18T13:41:03.977528Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"try:\n    rna_enrich_cl1 = gp.enrichr(top100_rna_non_overlap_clusters.loc[top100_rna_non_overlap_clusters.cluster == 0].index.tolist(), \n                            gene_sets=gene_sets, organism=\"Human\", cutoff=0.01)\nexcept:\n    rna_enrich_cl1 = gp.enrichr(top100_rna_non_overlap_clusters.loc[top100_rna_non_overlap_clusters.cluster == 0].index.tolist(), \n                            gene_sets=gene_sets, organism=\"Human\", cutoff=0.01)","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:41:05.877059Z","iopub.execute_input":"2023-04-18T13:41:05.878413Z","iopub.status.idle":"2023-04-18T13:41:11.063752Z","shell.execute_reply.started":"2023-04-18T13:41:05.878351Z","shell.execute_reply":"2023-04-18T13:41:11.062212Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"display(rna_enrich_cl1.results.sort_values(\"Adjusted P-value\").head(10))","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:41:12.262966Z","iopub.execute_input":"2023-04-18T13:41:12.263855Z","iopub.status.idle":"2023-04-18T13:41:12.283809Z","shell.execute_reply.started":"2023-04-18T13:41:12.263805Z","shell.execute_reply":"2023-04-18T13:41:12.282512Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"try:\n    rna_enrich_cl2 = gp.enrichr(top100_rna_non_overlap_clusters.loc[top100_rna_non_overlap_clusters.cluster == 1].index.tolist(), \n                            gene_sets=gene_sets, organism=\"Human\", cutoff=0.01)\nexcept:\n    rna_enrich_cl2 = gp.enrichr(top100_rna_non_overlap_clusters.loc[top100_rna_non_overlap_clusters.cluster == 1].index.tolist(), \n                            gene_sets=gene_sets, organism=\"Human\", cutoff=0.01)\n","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:41:13.470348Z","iopub.execute_input":"2023-04-18T13:41:13.471629Z","iopub.status.idle":"2023-04-18T13:41:16.962881Z","shell.execute_reply.started":"2023-04-18T13:41:13.471579Z","shell.execute_reply":"2023-04-18T13:41:16.961555Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"display(rna_enrich_cl2.results.sort_values(\"Adjusted P-value\").head(10))","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:41:21.054662Z","iopub.execute_input":"2023-04-18T13:41:21.055255Z","iopub.status.idle":"2023-04-18T13:41:21.081136Z","shell.execute_reply.started":"2023-04-18T13:41:21.055205Z","shell.execute_reply":"2023-04-18T13:41:21.079710Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Clustering of non-overlap proteins","metadata":{}},{"cell_type":"code","source":"sns.clustermap(df_rna[top100_prot_non_overlap].corr()) ","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:41:25.380769Z","iopub.execute_input":"2023-04-18T13:41:25.381227Z","iopub.status.idle":"2023-04-18T13:41:26.264038Z","shell.execute_reply.started":"2023-04-18T13:41:25.381178Z","shell.execute_reply":"2023-04-18T13:41:26.262449Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"clustering_prot = AgglomerativeClustering(linkage=\"average\", n_clusters=2, \n                                     compute_full_tree=True)\nclustering_prot.fit_predict(df_rna[top100_prot_non_overlap].corr())\ntop100_prot_non_overlap_clusters = get_cluster_table(df_rna[top100_prot_non_overlap].corr().index, clustering_prot)","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:41:27.522711Z","iopub.execute_input":"2023-04-18T13:41:27.523193Z","iopub.status.idle":"2023-04-18T13:41:27.714378Z","shell.execute_reply.started":"2023-04-18T13:41:27.523149Z","shell.execute_reply":"2023-04-18T13:41:27.713283Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"try:\n    prot_enrich_cl1 = gp.enrichr(top100_prot_non_overlap_clusters.loc[top100_prot_non_overlap_clusters.cluster == 0].index.tolist(), \n                            gene_sets=gene_sets, organism=\"Human\", cutoff=0.01)\nexcept:\n    prot_enrich_cl1 = gp.enrichr(top100_prot_non_overlap_clusters.loc[top100_prot_non_overlap_clusters.cluster == 0].index.tolist(), \n                            gene_sets=gene_sets, organism=\"Human\", cutoff=0.01)","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:41:30.177790Z","iopub.execute_input":"2023-04-18T13:41:30.178913Z","iopub.status.idle":"2023-04-18T13:41:35.909683Z","shell.execute_reply.started":"2023-04-18T13:41:30.178863Z","shell.execute_reply":"2023-04-18T13:41:35.908685Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"display(prot_enrich_cl1.results.sort_values(\"Adjusted P-value\").head(10))","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:41:35.912315Z","iopub.execute_input":"2023-04-18T13:41:35.912638Z","iopub.status.idle":"2023-04-18T13:41:35.932960Z","shell.execute_reply.started":"2023-04-18T13:41:35.912608Z","shell.execute_reply":"2023-04-18T13:41:35.931621Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"try:\n    prot_enrich_cl2 = gp.enrichr(top100_prot_non_overlap_clusters.loc[top100_prot_non_overlap_clusters.cluster == 1].index.tolist(), \n                            gene_sets=gene_sets, organism=\"Human\", cutoff=0.01)\nexcept:\n    prot_enrich_cl2 = gp.enrichr(top100_prot_non_overlap_clusters.loc[top100_prot_non_overlap_clusters.cluster == 1].index.tolist(), \n                            gene_sets=gene_sets, organism=\"Human\", cutoff=0.01)","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:41:35.934811Z","iopub.execute_input":"2023-04-18T13:41:35.935303Z","iopub.status.idle":"2023-04-18T13:41:39.494713Z","shell.execute_reply.started":"2023-04-18T13:41:35.935255Z","shell.execute_reply":"2023-04-18T13:41:39.493765Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"display(prot_enrich_cl2.results.sort_values(\"Adjusted P-value\").head(10))","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:41:39.497753Z","iopub.execute_input":"2023-04-18T13:41:39.501420Z","iopub.status.idle":"2023-04-18T13:41:39.526358Z","shell.execute_reply.started":"2023-04-18T13:41:39.501345Z","shell.execute_reply":"2023-04-18T13:41:39.524728Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# RNA vs. Protein correlation scatter plot","metadata":{}},{"cell_type":"code","source":"df_importances_spearman_abs = df_importances_spearman.abs()\ndf_importances_spearman_abs[\"label\"] = None\ndf_importances_spearman_abs.loc[df_importances_spearman_abs[\"CD32 Prot Spearman Corr\"] > 0.1, \"label\"] = \"protein\"\ndf_importances_spearman_abs.loc[df_importances_spearman_abs[\"FCGR2A RNA Spearman Corr\"] > 0.1, \"label\"] = \"rna\"\ndf_importances_spearman_abs.loc[(df_importances_spearman_abs[\"FCGR2A RNA Spearman Corr\"] > 0.1) & \n                                (df_importances_spearman_abs[\"CD32 Prot Spearman Corr\"] > 0.1), \"label\"] = \"both\"\n        ","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:41:44.768110Z","iopub.execute_input":"2023-04-18T13:41:44.768975Z","iopub.status.idle":"2023-04-18T13:41:44.782136Z","shell.execute_reply.started":"2023-04-18T13:41:44.768931Z","shell.execute_reply":"2023-04-18T13:41:44.780887Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"fig = plt.figure(figsize = ( 10,10 ))\nsns.scatterplot(df_importances_spearman_abs.loc[df_importances_spearman_abs.label.isna() == False], \n                 x=\"CD32 Prot Spearman Corr\", y=\"FCGR2A RNA Spearman Corr\", hue=\"label\")\nplt.xlim([0.0, 0.3])\nplt.ylim([0.0, 0.3])","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:41:45.511664Z","iopub.execute_input":"2023-04-18T13:41:45.512864Z","iopub.status.idle":"2023-04-18T13:41:46.064740Z","shell.execute_reply.started":"2023-04-18T13:41:45.512818Z","shell.execute_reply":"2023-04-18T13:41:46.063402Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Non-overlap RNA enrichmment","metadata":{}},{"cell_type":"code","source":"df_importances_spearman[df_importances_spearman_abs.label == \"rna\"].sort_values(\"FCGR2A RNA Spearman Corr\", ascending=False).head(20)","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:43:26.454245Z","iopub.execute_input":"2023-04-18T13:43:26.455693Z","iopub.status.idle":"2023-04-18T13:43:26.483270Z","shell.execute_reply.started":"2023-04-18T13:43:26.455607Z","shell.execute_reply":"2023-04-18T13:43:26.481713Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"try:\n    enrich = gp.enrichr(df_importances_spearman_abs.loc[df_importances_spearman_abs.label == \"rna\"].index.tolist(), \n                            gene_sets=gene_sets, organism=\"Human\", cutoff=0.01)\nexcept:\n    enrich = gp.enrichr(df_importances_spearman_abs.loc[df_importances_spearman_abs.label == \"rna\"].index.tolist(), \n                            gene_sets=gene_sets, organism=\"Human\", cutoff=0.01)","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:43:38.533578Z","iopub.execute_input":"2023-04-18T13:43:38.534050Z","iopub.status.idle":"2023-04-18T13:43:43.589502Z","shell.execute_reply.started":"2023-04-18T13:43:38.534008Z","shell.execute_reply":"2023-04-18T13:43:43.588102Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"display(enrich.results.sort_values(\"Adjusted P-value\").head(10))","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:43:43.592281Z","iopub.execute_input":"2023-04-18T13:43:43.592802Z","iopub.status.idle":"2023-04-18T13:43:43.616117Z","shell.execute_reply.started":"2023-04-18T13:43:43.592747Z","shell.execute_reply":"2023-04-18T13:43:43.614713Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Non-overlap protein enrichmment","metadata":{}},{"cell_type":"code","source":"df_importances_spearman[df_importances_spearman_abs.label == \"protein\"].sort_values(\"CD32 Prot Spearman Corr\", ascending=False).head(20)","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:43:57.220809Z","iopub.execute_input":"2023-04-18T13:43:57.221231Z","iopub.status.idle":"2023-04-18T13:43:57.241394Z","shell.execute_reply.started":"2023-04-18T13:43:57.221197Z","shell.execute_reply":"2023-04-18T13:43:57.240051Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"try:\n    enrich = gp.enrichr(df_importances_spearman_abs.loc[df_importances_spearman_abs.label == \"protein\"].index.tolist(), \n                            gene_sets=gene_sets, organism=\"Human\", cutoff=0.01)\nexcept:\n    enrich = gp.enrichr(df_importances_spearman_abs.loc[df_importances_spearman_abs.label == \"protein\"].index.tolist(), \n                            gene_sets=gene_sets, organism=\"Human\", cutoff=0.01)","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:47:06.726287Z","iopub.execute_input":"2023-04-18T13:47:06.726767Z","iopub.status.idle":"2023-04-18T13:47:13.690056Z","shell.execute_reply.started":"2023-04-18T13:47:06.726728Z","shell.execute_reply":"2023-04-18T13:47:13.688932Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"display(enrich.results.sort_values(\"Adjusted P-value\").head(10))","metadata":{"execution":{"iopub.status.busy":"2023-04-18T13:47:13.691926Z","iopub.execute_input":"2023-04-18T13:47:13.692290Z","iopub.status.idle":"2023-04-18T13:47:13.715020Z","shell.execute_reply.started":"2023-04-18T13:47:13.692256Z","shell.execute_reply":"2023-04-18T13:47:13.713569Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}