{"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\nGrandmaster Silogram kindly shared lists of his important features for the past competition on single cell data.\nhttps://www.kaggle.com/competitions/open-problems-multimodal/discussion/366455\n\nWe can see quite some biology from them - e.g. -  CD36 protein is higly activated in erythroid like cells - similar to red blood cells - so you can see genes like HBD, HBB, HBA1 as important features - that various forms  of the hemoglobin  - so quite as expected from biology.\n\nAlso look on the so-called genes enrichment analysis with KEGG pathways - again biologically reasonable results,\nand compared with importances by other methods - they are quite consistent.\n\nThat is a first look - we need some time to get more insights from the data. \n\n","metadata":{}},{"cell_type":"markdown","source":"# Load and look on data","metadata":{}},{"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 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        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","execution":{"iopub.status.busy":"2023-01-09T21:31:41.796778Z","iopub.execute_input":"2023-01-09T21:31:41.797337Z","iopub.status.idle":"2023-01-09T21:31:41.886139Z","shell.execute_reply.started":"2023-01-09T21:31:41.797218Z","shell.execute_reply":"2023-01-09T21:31:41.884742Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\nfn2 = '/kaggle/input/featurespickle/CSVfeaturesCSV.csv'\nfn = '/kaggle/input/citeq-feature-importance/features_04.10.csv'\ndf = pd.read_csv(fn)\ndf","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:10:06.858165Z","iopub.execute_input":"2023-01-09T20:10:06.858514Z","iopub.status.idle":"2023-01-09T20:10:07.172385Z","shell.execute_reply.started":"2023-01-09T20:10:06.858483Z","shell.execute_reply":"2023-01-09T20:10:07.170948Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Look on feature importances for CD36 ","metadata":{}},{"cell_type":"code","source":"df[df['target'] == 'CD36']","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:10:07.174462Z","iopub.execute_input":"2023-01-09T20:10:07.175018Z","iopub.status.idle":"2023-01-09T20:10:07.218144Z","shell.execute_reply.started":"2023-01-09T20:10:07.174946Z","shell.execute_reply":"2023-01-09T20:10:07.217074Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df[df['target'] == 'CD36'].to_csv('Silogram CD36 Feature importances.csv')","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:10:07.221148Z","iopub.execute_input":"2023-01-09T20:10:07.222575Z","iopub.status.idle":"2023-01-09T20:10:07.252251Z","shell.execute_reply.started":"2023-01-09T20:10:07.222520Z","shell.execute_reply":"2023-01-09T20:10:07.251183Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df[df['target'] == 'CD36'].head(50)","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:10:07.253736Z","iopub.execute_input":"2023-01-09T20:10:07.255010Z","iopub.status.idle":"2023-01-09T20:10:07.293059Z","shell.execute_reply.started":"2023-01-09T20:10:07.254929Z","shell.execute_reply":"2023-01-09T20:10:07.291731Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import matplotlib.pyplot as plt\nv = df[df['target'] == 'CD36'].iloc[:50]['Importance'].values\n#print(v)\nfig = plt.figure(figsize = (20,5))\nplt.plot(v , '*-' )\nplt.show()\nfig = plt.figure(figsize = (20,5))\nplt.bar( x= range(len(v)), height =  v)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:10:07.294813Z","iopub.execute_input":"2023-01-09T20:10:07.296157Z","iopub.status.idle":"2023-01-09T20:10:07.856087Z","shell.execute_reply.started":"2023-01-09T20:10:07.296105Z","shell.execute_reply":"2023-01-09T20:10:07.854712Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import matplotlib.pyplot as plt\nv = df[df['target'] == 'CD36'].iloc[1:50]['Importance'].values\n#print(v)\nfig = plt.figure(figsize = (20,5))\nplt.plot(v , '*-' )\nplt.show()\nfig = plt.figure(figsize = (20,5))\nplt.bar( x= range(len(v)), height =  v)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:10:07.857825Z","iopub.execute_input":"2023-01-09T20:10:07.858354Z","iopub.status.idle":"2023-01-09T20:10:08.383738Z","shell.execute_reply.started":"2023-01-09T20:10:07.858306Z","shell.execute_reply":"2023-01-09T20:10:08.382494Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import matplotlib.pyplot as plt\nv = df[df['target'] == 'CD36'].iloc[10:100]['Importance'].values\n#print(v)\nfig = plt.figure(figsize = (20,5))\nplt.plot(v , '*-' )\nplt.show()\nfig = plt.figure(figsize = (20,5))\nplt.bar( x= range(len(v)), height =  v)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:10:08.385608Z","iopub.execute_input":"2023-01-09T20:10:08.386457Z","iopub.status.idle":"2023-01-09T20:10:09.230069Z","shell.execute_reply.started":"2023-01-09T20:10:08.386409Z","shell.execute_reply":"2023-01-09T20:10:09.228565Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import matplotlib.pyplot as plt\nv = df[df['target'] == 'CD36'].iloc[10:1000]['Importance'].values\n#print(v)\nv = np.log(1+v)\nfig = plt.figure(figsize = (20,5))\nplt.plot(v , '*-' )\nplt.show()\nfig = plt.figure(figsize = (20,5))\nplt.bar( x= range(len(v)), height =  v)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:10:09.231731Z","iopub.execute_input":"2023-01-09T20:10:09.232238Z","iopub.status.idle":"2023-01-09T20:10:11.584903Z","shell.execute_reply.started":"2023-01-09T20:10:09.232188Z","shell.execute_reply":"2023-01-09T20:10:11.583647Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Genes sets enrichment analysis - what pathways are related to important genes ","metadata":{}},{"cell_type":"code","source":"!pip install gseapy","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:10:11.588833Z","iopub.execute_input":"2023-01-09T20:10:11.589366Z","iopub.status.idle":"2023-01-09T20:10:25.819913Z","shell.execute_reply.started":"2023-01-09T20:10:11.589329Z","shell.execute_reply":"2023-01-09T20:10:25.818312Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import gseapy as gp","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:10:25.821841Z","iopub.execute_input":"2023-01-09T20:10:25.822264Z","iopub.status.idle":"2023-01-09T20:10:26.201700Z","shell.execute_reply.started":"2023-01-09T20:10:25.822228Z","shell.execute_reply":"2023-01-09T20:10:26.200370Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## 100 genes ","metadata":{}},{"cell_type":"code","source":"list_genes = df[df['target'] == 'CD36']['Feature'].tolist()\nlist_genes = list_genes[:100]\nlist_genes = [t.split('_')[1] for t in list_genes]\nprint(list_genes[:100])","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:10:26.205629Z","iopub.execute_input":"2023-01-09T20:10:26.206173Z","iopub.status.idle":"2023-01-09T20:10:26.232833Z","shell.execute_reply.started":"2023-01-09T20:10:26.206125Z","shell.execute_reply":"2023-01-09T20:10:26.231079Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\nenr = gp.enrichr(\n    gene_list=list_genes,\n    gene_sets=['KEGG_2016','KEGG_2021_Human'],\n    organism='human',\n    outdir=None,\n)","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:10:26.235339Z","iopub.execute_input":"2023-01-09T20:10:26.236450Z","iopub.status.idle":"2023-01-09T20:10:29.503952Z","shell.execute_reply.started":"2023-01-09T20:10:26.236411Z","shell.execute_reply":"2023-01-09T20:10:29.502604Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"enr.results.head(50)","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:10:29.505509Z","iopub.execute_input":"2023-01-09T20:10:29.505873Z","iopub.status.idle":"2023-01-09T20:10:29.536681Z","shell.execute_reply.started":"2023-01-09T20:10:29.505842Z","shell.execute_reply":"2023-01-09T20:10:29.535666Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## 500 genes ","metadata":{}},{"cell_type":"code","source":"%%time\nlist_genes = df[df['target'] == 'CD36']['Feature'].tolist()\nlist_genes = list_genes[:500]\nlist_genes = [t.split('_')[1] for t in list_genes]\nprint(list_genes[:100])\n\nenr = gp.enrichr(\n    gene_list=list_genes,\n    gene_sets=['KEGG_2016','KEGG_2021_Human'],\n    organism='human',\n    outdir=None,\n)\n\nenr.results.head(50)","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:10:29.538077Z","iopub.execute_input":"2023-01-09T20:10:29.539323Z","iopub.status.idle":"2023-01-09T20:10:33.799120Z","shell.execute_reply.started":"2023-01-09T20:10:29.539274Z","shell.execute_reply":"2023-01-09T20:10:33.797823Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Enrichment for transcription factors ","metadata":{}},{"cell_type":"code","source":"list_tf_related_signatures = [ \n#'ChEA_2013',  'ChEA_2015', 'ChEA_2016', \n'ChEA_2022',  \n'ENCODE_TF_ChIP-seq_2015', 'ENCODE_and_ChEA_Consensus_TFs_from_ChIP-X', \n'TF-LOF_Expression_from_GEO', 'TF_Perturbations_Followed_by_Expression', \n'TRANSFAC_and_JASPAR_PWMs',\n'TRRUST_Transcription_Factors_2019', \n'ARCHS4_TFs_Coexp', 'Enrichr_Submissions_TF-Gene_Coocurrence',\n'Transcription_Factor_PPIs',\n] ","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:10:33.800922Z","iopub.execute_input":"2023-01-09T20:10:33.801312Z","iopub.status.idle":"2023-01-09T20:10:33.806914Z","shell.execute_reply.started":"2023-01-09T20:10:33.801277Z","shell.execute_reply":"2023-01-09T20:10:33.805576Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Enrichr_Submissions_TF-Gene_Coocurrence\nprint(len(list_genes), list_genes[:50])\nfor sig in list_tf_related_signatures: #  ['ChEA_2022',     'ENCODE_TF_ChIP-seq_2015', 'ENCODE_and_ChEA_Consensus_TFs_from_ChIP-X', \n#     'TF-LOF_Expression_from_GEO', 'TF_Perturbations_Followed_by_Expression',     'TRANSFAC_and_JASPAR_PWMs',\n#     'TRRUST_Transcription_Factors_2019',   ]:\n    \n    print(); print(sig)\n    if 'ChEA' in sig:\n        print(' \"ChEA\" - see: https://maayanlab.cloud/Harmonizome/dataset/CHEA+Transcription+Factor+Targets ' )\n        \n    enr = gp.enrichr(\n        gene_list=list_genes, # 'MSigDB_Oncogenic_Signatures',\n        gene_sets= sig,\n        organism='human',\n        outdir=None,\n    )\n    col = 'Adjusted P-value' #P-value # Adjusted P-value - with Bonferoni type correction for checking multiple lists, while just \"P-value\" - is without\n    m = enr.results[col] < 0.05\n    print(m.sum(), 'out of:', len( enr.results ) )\n    display( enr.results[m].sort_values(col).head(20) )\n#     display( enr.results.sort_values('Adjusted P-value').head(20) )","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:10:33.808925Z","iopub.execute_input":"2023-01-09T20:10:33.809899Z","iopub.status.idle":"2023-01-09T20:10:58.081102Z","shell.execute_reply.started":"2023-01-09T20:10:33.809858Z","shell.execute_reply":"2023-01-09T20:10:58.079801Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Compare with the other feature importances - correlations , catboost, etc..","metadata":{}},{"cell_type":"code","source":"fn = '/kaggle/input/research-project-01-around-multimodal-singlecell/CD36importances.csv'\ndf2 = pd.read_csv(fn)\ndf2","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:10:58.082899Z","iopub.execute_input":"2023-01-09T20:10:58.083447Z","iopub.status.idle":"2023-01-09T20:10:58.116791Z","shell.execute_reply.started":"2023-01-09T20:10:58.083402Z","shell.execute_reply":"2023-01-09T20:10:58.115616Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df2.columns = list( df2.iloc[0,:])\ndf2=df2.iloc[1:,:]\ndf2","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:10:58.118446Z","iopub.execute_input":"2023-01-09T20:10:58.119111Z","iopub.status.idle":"2023-01-09T20:10:58.142021Z","shell.execute_reply.started":"2023-01-09T20:10:58.119074Z","shell.execute_reply":"2023-01-09T20:10:58.140798Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_stat = pd.DataFrame(); IX = 0\n\nfor N in [100,200,500,1000]:\n    list_genes = df[df['target'] == 'CD36']['Feature'].tolist()\n    list_genes = list_genes[:N]\n    \n\n    print('Interesection of Silogram  with others methods')\n    for col in df2.columns:\n        s = set([t.split('_')[1] for t in list_genes]) & set( [str(t).split('_')[-1] for t in df2[col]] )\n        #print(col, len(s))\n        df_stat.loc[col, 'Intersection with Silogram '+str(N)] = len(s)\n        df_stat.loc[col, 'Out of'] = len(set( df2[col]) )\n\ndf_stat    ","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:10:58.143076Z","iopub.execute_input":"2023-01-09T20:10:58.143427Z","iopub.status.idle":"2023-01-09T20:10:58.264904Z","shell.execute_reply.started":"2023-01-09T20:10:58.143396Z","shell.execute_reply":"2023-01-09T20:10:58.263917Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Analysis for other proteins ","metadata":{}},{"cell_type":"markdown","source":"## Mouse IgG2a","metadata":{}},{"cell_type":"code","source":"fn = '/kaggle/input/silogram-fi-by-cd/Mouse-IgG2a_protein_feature_importances.csv'\ndf = pd.read_csv(fn,index_col = 0)\ndf","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:10:58.266242Z","iopub.execute_input":"2023-01-09T20:10:58.266546Z","iopub.status.idle":"2023-01-09T20:10:58.288821Z","shell.execute_reply.started":"2023-01-09T20:10:58.266519Z","shell.execute_reply":"2023-01-09T20:10:58.287574Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"l = []\nfor k in range(df.shape[0]):\n    s = set( df.iloc[:k,0]) & set( df.iloc[:k,1])\n    l.append(len(s))    ","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:10:58.290154Z","iopub.execute_input":"2023-01-09T20:10:58.290589Z","iopub.status.idle":"2023-01-09T20:10:58.320515Z","shell.execute_reply.started":"2023-01-09T20:10:58.290512Z","shell.execute_reply":"2023-01-09T20:10:58.319629Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import matplotlib.pyplot as plt\nplt.plot(l)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:10:58.321788Z","iopub.execute_input":"2023-01-09T20:10:58.322417Z","iopub.status.idle":"2023-01-09T20:10:58.536892Z","shell.execute_reply.started":"2023-01-09T20:10:58.322380Z","shell.execute_reply":"2023-01-09T20:10:58.535735Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_stat = pd.DataFrame()\nll = os.listdir('/kaggle/input/silogram-fi-by-cd/' )\nc = 0\nfor f in ll:\n    if 'csv' not in f: continue\n    c += 1\n    #print(f)\n    fn = '/kaggle/input/silogram-fi-by-cd/' + f\n    #fn = '/kaggle/input/silogram-fi-by-cd/Mouse-IgG2a_protein_feature_importances.csv'\n    df = pd.read_csv(fn,index_col = 0)\n    \n    cd_name = f.split('_')[0]\n    k = df.shape[0]\n    s = set( df.iloc[:k,0]) & set( df.iloc[:k,1])\n    df_stat.loc[cd_name,'Intersection Silogram 1,2'] = len(s)\n    \ncol = df_stat.columns[0]\ndf_stat = df_stat.sort_values(col, ascending = False)\ndisplay(df_stat.head(10) )\ndisplay(df_stat.tail(10) )\n\n","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:10:58.538599Z","iopub.execute_input":"2023-01-09T20:10:58.540057Z","iopub.status.idle":"2023-01-09T20:10:59.454178Z","shell.execute_reply.started":"2023-01-09T20:10:58.540001Z","shell.execute_reply":"2023-01-09T20:10:59.452935Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"col = df_stat.columns[0]\nplt.hist(df_stat[col], bins = 30)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:10:59.455851Z","iopub.execute_input":"2023-01-09T20:10:59.456790Z","iopub.status.idle":"2023-01-09T20:10:59.693287Z","shell.execute_reply.started":"2023-01-09T20:10:59.456747Z","shell.execute_reply":"2023-01-09T20:10:59.692069Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#%%time\ndf_y = pd.read_hdf('/kaggle/input/open-problems-multimodal/train_cite_targets.h5')\ndisplay(df_y)","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:10:59.694826Z","iopub.execute_input":"2023-01-09T20:10:59.696319Z","iopub.status.idle":"2023-01-09T20:11:00.612443Z","shell.execute_reply.started":"2023-01-09T20:10:59.696267Z","shell.execute_reply":"2023-01-09T20:11:00.611041Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"s = set(df_y.columns) & set(df_stat.index)\nlen(s)\n","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:11:00.614228Z","iopub.execute_input":"2023-01-09T20:11:00.614604Z","iopub.status.idle":"2023-01-09T20:11:00.621633Z","shell.execute_reply.started":"2023-01-09T20:11:00.614573Z","shell.execute_reply":"2023-01-09T20:11:00.620578Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"t = df_y.mean(axis = 0)\nt.name = 'Protein mean'\ndf_stat2 = df_stat.join(t)\n\nt = df_y.std(axis = 0)\nt.name = 'Protein std'\ndf_stat2 = df_stat2.join(t)\n\nt = df_y.var(axis = 0)\nt.name = 'Protein var'\ndf_stat2 = df_stat2.join(t)\n\ndf_stat2","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:11:00.626934Z","iopub.execute_input":"2023-01-09T20:11:00.628057Z","iopub.status.idle":"2023-01-09T20:11:00.906424Z","shell.execute_reply.started":"2023-01-09T20:11:00.628014Z","shell.execute_reply":"2023-01-09T20:11:00.905092Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_stat2.to_csv('CD_proteins_stat.csv')","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:11:00.908044Z","iopub.execute_input":"2023-01-09T20:11:00.908425Z","iopub.status.idle":"2023-01-09T20:11:00.917325Z","shell.execute_reply.started":"2023-01-09T20:11:00.908392Z","shell.execute_reply":"2023-01-09T20:11:00.916046Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_stat2.corr()","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:11:00.918860Z","iopub.execute_input":"2023-01-09T20:11:00.921051Z","iopub.status.idle":"2023-01-09T20:11:00.936795Z","shell.execute_reply.started":"2023-01-09T20:11:00.920995Z","shell.execute_reply":"2023-01-09T20:11:00.935452Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Linear growth pattern between intersection of two top100 lists - indicitation of their non-randomness.\n\nFor random lists the growth would be quadractic - parabola x^2/N \nhttps://mathoverflow.net/q/437545/10446\n\nSee also examples here:\nhttps://www.kaggle.com/code/alexandervc/mmscel-correlations-cd-vs-rna-v2?scriptVersionId=115306085&cellId=90\n\n\nThe point where the growth deviates from the linear - is an indication how many features one can take. \n\n(It is just indication - not a competely \"proved\" scheme for the moment). ","metadata":{}},{"cell_type":"code","source":"ll = os.listdir('/kaggle/input/silogram-fi-by-cd/' )\n\nn_x_subplots = 6\n\nc = 0\nfor f in ll:\n    if 'csv' not in f: continue\n    \n    if c % n_x_subplots == 0:\n        if c > 0:\n            plt.show()\n        fig = plt.figure(figsize = (20,4) ); c = 0\n        #plt.suptitle(str(k),fontsize = 20 )# str_data_inf + ' n_cells: ' +str(mask.sum()) + ' RED > median expression, BLUE <= median ' )# +' ' + cell_type +' ' + drug )\n    c += 1; fig.add_subplot(1,n_x_subplots ,c)\n    \n    print(f)\n    fn = '/kaggle/input/silogram-fi-by-cd/' + f\n    #fn = '/kaggle/input/silogram-fi-by-cd/Mouse-IgG2a_protein_feature_importances.csv'\n    df = pd.read_csv(fn,index_col = 0)\n    \n    l = []\n    for k in range(df.shape[0]):\n        s = set( df.iloc[:k,0]) & set( df.iloc[:k,1])\n        l.append(len(s))    \n\n    plt.plot(l)\n    plt.title(f.split('_')[0] )\n    plt.grid()\n    \n    if c >= 10: break\n        \nplt.show()        \n","metadata":{"execution":{"iopub.status.busy":"2023-01-09T20:11:00.938407Z","iopub.execute_input":"2023-01-09T20:11:00.938848Z","iopub.status.idle":"2023-01-09T20:11:22.202638Z","shell.execute_reply.started":"2023-01-09T20:11:00.938813Z","shell.execute_reply":"2023-01-09T20:11:22.201493Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\nfn2 = '/kaggle/input/featurespickle/CSVfeaturesCSV.csv'\nfn = '/kaggle/input/citeq-feature-importance/features_04.10.csv'\ndf = pd.read_csv(fn)\ndisplay(df)\ndf2 = pd.read_csv(fn2, index_col = 0)\ndisplay(df2)\n","metadata":{"execution":{"iopub.status.busy":"2023-01-09T21:38:38.197207Z","iopub.execute_input":"2023-01-09T21:38:38.197665Z","iopub.status.idle":"2023-01-09T21:38:42.864134Z","shell.execute_reply.started":"2023-01-09T21:38:38.197603Z","shell.execute_reply":"2023-01-09T21:38:42.862895Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"m = df['target'] == 'CD36'\nprint( m.sum() )\nm2 = df2['target'] == 'CD36'\nprint( m2.sum() )\ndf_importances = df2[m2].set_index('Feature').drop('target', axis = 1)\ndf_importances.columns = ['Trial0']\ndf_importances = df_importances.join( df[m].set_index('Feature').drop('target', axis = 1), how = 'inner' ) \ndf_importances.columns = ['Trial0','Trial1']\nprint( df_importances['Trial0'].isnull().sum() )\nprint( df_importances['Trial1'].isnull().sum() )\ndf_importances","metadata":{"execution":{"iopub.status.busy":"2023-01-09T21:43:29.431754Z","iopub.execute_input":"2023-01-09T21:43:29.432145Z","iopub.status.idle":"2023-01-09T21:43:29.673696Z","shell.execute_reply.started":"2023-01-09T21:43:29.432114Z","shell.execute_reply":"2023-01-09T21:43:29.672460Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"display(df_importances.describe())\ndisplay(df_importances.corr())\ndisplay( (-df_importances).rank().head(20) )\n","metadata":{"execution":{"iopub.status.busy":"2023-01-09T21:47:37.345234Z","iopub.execute_input":"2023-01-09T21:47:37.345692Z","iopub.status.idle":"2023-01-09T21:47:37.389082Z","shell.execute_reply.started":"2023-01-09T21:47:37.345622Z","shell.execute_reply":"2023-01-09T21:47:37.387764Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import matplotlib.pyplot as plt\nimport seaborn as sns\n","metadata":{"execution":{"iopub.status.busy":"2023-01-09T21:46:05.182049Z","iopub.execute_input":"2023-01-09T21:46:05.182474Z","iopub.status.idle":"2023-01-09T21:46:05.764226Z","shell.execute_reply.started":"2023-01-09T21:46:05.182440Z","shell.execute_reply":"2023-01-09T21:46:05.763008Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time \ndict_counts = {}\n\ndict_sorted = {} # pd.DataFrame()\nfor col in df_importances.columns:\n    dict_sorted[col] = df_importances[col].abs().sort_values(ascending = False )\n\nfor k in range(0, df_importances.shape[0]):\n    if (k%5000 == 1): print(k)\n    col = df_importances.columns[0]\n    s = set(  dict_sorted[col].index[:k] )\n    for col in df_importances.columns:\n        \n        s = s & set(dict_sorted[col].index[:k] )\n        if col not in dict_counts.keys(): \n            dict_counts[col] = []\n        dict_counts[col].append(len(s))\n\nplt.figure(figsize = (20,10))\nfor i,col in enumerate(dict_counts):\n    if len(dict_counts ) >= 50:\n        if (i%10) != 0: continue\n    plt.plot(dict_counts[col], label = 'intersection till: ' + col )\nplt.legend(fontsize = 20 )\nplt.grid()\nplt.show()\n\nplt.figure(figsize = (20,10))\nfor i,col in enumerate(dict_counts):\n    if len(dict_counts ) >= 50:\n        if (i%10) != 0: continue\n    plt.plot(dict_counts[col], label = 'intersection till: ' + col )\nplt.legend(fontsize = 20 )\nplt.xlim([0,200])\nplt.grid()\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-01-09T21:49:26.442983Z","iopub.execute_input":"2023-01-09T21:49:26.444393Z","iopub.status.idle":"2023-01-09T21:49:27.295677Z","shell.execute_reply.started":"2023-01-09T21:49:26.444350Z","shell.execute_reply":"2023-01-09T21:49:27.294510Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Plots of parabolic approximation for full list of features","metadata":{}},{"cell_type":"code","source":"for i,col in  enumerate(dict_counts): # = '3 32606 EryP'\n    if i == 0: continue \n    if len(dict_counts ) >= 50:\n        if (i%10) != 0: continue \n\n    ll = dict_counts[col]\n    plt.figure(figsize = (20,6))\n    plt.plot(ll, label = 'intersection till: ' + col)\n    y = np.array(ll)\n    x = np.arange(len(y))\n    p = np.polyfit(x,y,2)\n    print(p)\n    plt.plot(np.polyval(p,x) , label = 'quadratic approx'  )\n    plt.legend(fontsize = 20 )\n    plt.grid()\n    plt.show()","metadata":{"execution":{"iopub.status.busy":"2023-01-09T21:56:00.579067Z","iopub.execute_input":"2023-01-09T21:56:00.579600Z","iopub.status.idle":"2023-01-09T21:56:00.833775Z","shell.execute_reply.started":"2023-01-09T21:56:00.579560Z","shell.execute_reply":"2023-01-09T21:56:00.832466Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"list_slopes = []\nfor i, col in  enumerate(dict_counts): # = '3 32606 EryP'\n    if i == 0: continue \n    ll = dict_counts[col]\n    \n    y = np.array(ll)\n    x = np.arange(len(y))\n    \n    M = 100\n    p = np.polyfit(x[:M],y[:M],1)\n    print(p)\n    list_slopes.append(p[0])\n    \n    if len(dict_counts ) >= 50:\n        if (i%10) != 0: continue \n        \n    plt.figure(figsize = (20,6))\n    plt.plot(ll, label = 'intersection till: ' + col)\n    y = np.array(ll)\n    x = np.arange(len(y))\n    p = np.polyfit(x,y,2)\n    print(p)\n    plt.plot(np.polyval(p,x), label = 'quadratic approx'  )\n\n    M = 100\n    p = np.polyfit(x[:M],y[:M],1)\n    print(p)\n    plt.plot(np.polyval(p,x), label = 'linear approx' )\n\n    \n    if i < 3:\n        plt.ylim([0,300])\n    else:\n        plt.ylim([0,300])\n    plt.xlim([0,300])\n    plt.legend(fontsize = 20 )\n    plt.grid()\n    plt.show()","metadata":{"execution":{"iopub.status.busy":"2023-01-09T21:53:38.562071Z","iopub.execute_input":"2023-01-09T21:53:38.562695Z","iopub.status.idle":"2023-01-09T21:53:38.866511Z","shell.execute_reply.started":"2023-01-09T21:53:38.562643Z","shell.execute_reply":"2023-01-09T21:53:38.865064Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}