{"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\nFor all proteins we compare several simple Ridge based models for predicting it:\n\nUse only one feature - own RNA for protein e.g. for CD36 protein  we use CD36 rna - and predict by Ridge only based on that model\n\nUse cell type (one-hot encoded for prediction) \n\nUse cell type + own RNA\n\nUse top 100  PCA feautures as predictors \n\n\n#### Versions\n\n##### 1 - no allRNA features because to slow\n##### 3 - catboost on allRNA features, GPU\n##### 5 - SHAP catboost on important allRNA features (with PredictionValuesChange > 0.3), GPU\n(important features - with PredictionValuesChange > N) https://catboost.ai/en/docs/en/concepts/fstr#regular-feature-importance\n##### 6 - SHAP catboost on important allRNA features (with PredictionValuesChange > 0.1), GPU\n\n","metadata":{}},{"cell_type":"markdown","source":"# Key Params","metadata":{"execution":{"iopub.status.busy":"2022-12-18T21:18:08.728359Z","iopub.execute_input":"2022-12-18T21:18:08.729709Z","iopub.status.idle":"2022-12-18T21:18:08.735176Z","shell.execute_reply.started":"2022-12-18T21:18:08.729662Z","shell.execute_reply":"2022-12-18T21:18:08.733772Z"}}},{"cell_type":"code","source":"# target_name = 'CD36'\nn_splits_for_cross_valdition = 2 # if bigger, then for models on all rna (22K features) it may cause RAM crash and would be slower - \n\nrandom_state_cross_validation = 0 # random to be used in generating the folds - changing it we can study how stable are our results","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2023-01-11T15:38:56.356907Z","iopub.execute_input":"2023-01-11T15:38:56.357537Z","iopub.status.idle":"2023-01-11T15:38:56.380609Z","shell.execute_reply.started":"2023-01-11T15:38:56.357460Z","shell.execute_reply":"2023-01-11T15:38:56.379805Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Preparations","metadata":{}},{"cell_type":"code","source":"!pip install -q scanpy\nimport scanpy as sc\nimport anndata","metadata":{"execution":{"iopub.status.busy":"2023-01-11T15:38:56.386044Z","iopub.execute_input":"2023-01-11T15:38:56.386301Z","iopub.status.idle":"2023-01-11T15:39:14.665390Z","shell.execute_reply.started":"2023-01-11T15:38:56.386277Z","shell.execute_reply":"2023-01-11T15:39:14.664355Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\nfrom scipy import stats\n\nfrom sklearn.preprocessing import OneHotEncoder\nfrom sklearn.model_selection import cross_val_score, cross_val_predict, KFold\nfrom sklearn.metrics import r2_score, mean_squared_error, mean_absolute_error\n\nfrom sklearn.linear_model import Ridge, RidgeCV\nimport lightgbm as lgbm\nfrom catboost import CatBoostRegressor, Pool\nimport catboost\nimport shap\n\nimport time\nt0start = time.time()\n\nimport os, gc\nfrom tqdm import tqdm\n","metadata":{"execution":{"iopub.status.busy":"2023-01-11T15:39:14.667342Z","iopub.execute_input":"2023-01-11T15:39:14.668264Z","iopub.status.idle":"2023-01-11T15:39:17.204985Z","shell.execute_reply.started":"2023-01-11T15:39:14.668226Z","shell.execute_reply":"2023-01-11T15:39:17.204043Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Main results will be stored here: \ndf_scores = pd.DataFrame()","metadata":{"execution":{"iopub.status.busy":"2023-01-11T15:39:17.206593Z","iopub.execute_input":"2023-01-11T15:39:17.207373Z","iopub.status.idle":"2023-01-11T15:39:17.214001Z","shell.execute_reply.started":"2023-01-11T15:39:17.207335Z","shell.execute_reply":"2023-01-11T15:39:17.213013Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Load data for target","metadata":{}},{"cell_type":"code","source":"%%time\nfn = '/kaggle/input/citeseqscrnaseqproteins-challenge-neurips2021/GSE194122_openproblems_neurips2021_cite_BMMC_processed.h5ad'\nadata = sc.read(fn)\nprint(adata)\ncite_prots = adata[:, adata.var.feature_types == 'ADT']\ncite_genes = adata[:, adata.var.feature_types == 'GEX']\n# del adata\n# gc.collect()\n\ndf_y = pd.read_hdf('/kaggle/input/open-problems-multimodal/train_cite_targets.h5')\ndf_y","metadata":{"execution":{"iopub.status.busy":"2023-01-11T15:39:17.216636Z","iopub.execute_input":"2023-01-11T15:39:17.217245Z","iopub.status.idle":"2023-01-11T15:39:31.984627Z","shell.execute_reply.started":"2023-01-11T15:39:17.217201Z","shell.execute_reply":"2023-01-11T15:39:31.983591Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print(f'Number of intersected cite prots in 2021 and 2022: {len(set(cite_prots.var.index).intersection(set(df_y.columns)))}')\nprint(f'Number of prots in 2021: {cite_prots.var.index.shape[0]}')\nprint(f'Number of prots in 2022: {df_y.columns.shape[0]}')\nprint('Non intersected proteins from 2022:')\n[i for i in df_y.columns if i not in cite_prots.var.index]","metadata":{"execution":{"iopub.status.busy":"2023-01-11T15:39:31.986213Z","iopub.execute_input":"2023-01-11T15:39:31.986832Z","iopub.status.idle":"2023-01-11T15:39:32.000309Z","shell.execute_reply.started":"2023-01-11T15:39:31.986795Z","shell.execute_reply":"2023-01-11T15:39:31.999328Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_y = pd.DataFrame(cite_prots.X.A, cite_prots.obs.index, cite_prots.var.index)\n\ndf_scores.loc['Mean', df_y.columns]  = df_y.mean(axis = 0)\ndf_scores.loc['Std', df_y.columns]  = df_y.std(axis = 0)\ndel df_y\ngc.collect()","metadata":{"execution":{"iopub.status.busy":"2023-01-11T15:39:32.002635Z","iopub.execute_input":"2023-01-11T15:39:32.003047Z","iopub.status.idle":"2023-01-11T15:39:32.831724Z","shell.execute_reply.started":"2023-01-11T15:39:32.003000Z","shell.execute_reply":"2023-01-11T15:39:32.830606Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"#  Predicting based on -  own RNA ","metadata":{}},{"cell_type":"code","source":"def calc_lgb(\n    X, y, protein, df_scores, name_exp,\n    n_splits_for_cross_valdition = n_splits_for_cross_valdition, \n    random_state_cross_validation = random_state_cross_validation\n    ):\n    model = lgbm.LGBMRegressor(random_state = 0) # Ridge(alpha = alpha_selected )\n    kf = KFold(n_splits=n_splits_for_cross_valdition,  shuffle=True, random_state= random_state_cross_validation )\n    y_pred = cross_val_predict(model, X, y, cv=kf)\n    s = np.corrcoef(y,y_pred)[0,1]\n    df_scores.loc[f'Corr LGB {name_exp}', protein ]  = s\n    return df_scores\n\ndef calc_ridge(\n    X, y, protein, df_scores, name_exp,\n    ALPHAS = [1e-1, 1, 1e1,1e2,1e3,1e4,1e5],\n    n_splits_for_cross_valdition = n_splits_for_cross_valdition, \n    random_state_cross_validation = random_state_cross_validation\n    ):\n    \n    if ALPHAS:\n        model = RidgeCV(alphas=ALPHAS).fit(X, y)\n    else:\n        print('calc ridge with default params')\n        model = RidgeCV().fit(X, y)\n        \n    alpha_selected = model.alpha_\n    print('Optimal alpha found by cross-validation: ', alpha_selected )\n\n\n    model = Ridge(alpha = alpha_selected )\n    kf = KFold(n_splits=n_splits_for_cross_valdition,  shuffle=True, random_state= random_state_cross_validation )\n    y_pred = cross_val_predict(model, X, y, cv=kf)\n    s = r2_score(y,y_pred)\n    df_scores.loc[f'r2 Ridge {name_exp}', protein ]  = s\n    s = np.corrcoef(y,y_pred)[0,1]\n    df_scores.loc[f'Corr Ridge {name_exp}', protein ]  = s\n    s = mean_squared_error(y,y_pred, squared=False) # squared=False -> RMSE, not MSE \n    df_scores.loc[f'RMSE Ridge {name_exp}', protein ]  = s\n    return df_scores\n","metadata":{"execution":{"iopub.status.busy":"2023-01-11T15:39:32.833280Z","iopub.execute_input":"2023-01-11T15:39:32.833848Z","iopub.status.idle":"2023-01-11T15:39:32.847183Z","shell.execute_reply.started":"2023-01-11T15:39:32.833807Z","shell.execute_reply":"2023-01-11T15:39:32.845795Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\n# for protein in cite_prots.var.index:\n#     if protein in cite_genes.var.index:\n#         print(protein)\n#         X = pd.DataFrame(\n#             cite_genes[:, cite_genes.var.index == protein].X.A, \n#             cite_genes.obs.index, \n#             [protein]\n#         )\n#         y = pd.DataFrame(\n#             cite_prots[:, cite_prots.var.index == protein].X.A, \n#             cite_prots.obs.index, \n#             [protein]\n#         )\n        \n#         tmp_df = pd.DataFrame(\n#             {'rna': X.values.flatten(),\n#              'prot': y.values.flatten()}\n#         )\n#         # corr own RNA:\n#         df_scores.loc['Corr ownRNA', protein ]  = tmp_df.corr(method = 'pearson').values[0,1]\n#         df_scores.loc['Corr Spearman ownRNA', protein ] = tmp_df.corr(method = 'spearman').values[0,1]\n        \n#         df_scores = calc_lgb(X.values, y.values.flatten(), protein,  df_scores, 'ownRNA')\n#         df_scores = calc_ridge(X.values, y.values.flatten(), protein,  df_scores, 'ownRNA')\n        ","metadata":{"execution":{"iopub.status.busy":"2023-01-11T15:39:32.849079Z","iopub.execute_input":"2023-01-11T15:39:32.849442Z","iopub.status.idle":"2023-01-11T15:39:32.862310Z","shell.execute_reply.started":"2023-01-11T15:39:32.849397Z","shell.execute_reply":"2023-01-11T15:39:32.861102Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"#  Predicting based on - cell type information ","metadata":{}},{"cell_type":"code","source":"# %%time\n# df_meta = pd.read_csv('/kaggle/input/feature-shop-for-multimodal-singlecell-competition/_citeseq_meta_all_text_also.csv', index_col = 0 )\n# df_meta = df_meta[df_meta['Train0OrTest1']==0]\n# df_meta","metadata":{"execution":{"iopub.status.busy":"2023-01-11T15:39:32.864083Z","iopub.execute_input":"2023-01-11T15:39:32.864611Z","iopub.status.idle":"2023-01-11T15:39:32.874277Z","shell.execute_reply.started":"2023-01-11T15:39:32.864572Z","shell.execute_reply":"2023-01-11T15:39:32.873432Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# # X = df_meta[['HSC cell_type', 'EryP cell_type', 'NeuP cell_type','MasP cell_type', 'MkP cell_type', 'BP cell_type', 'MoP cell_type' ]]\n# # X \n# X = pd.get_dummies(adata.obs[['cell_type']])\n\n# for protein in cite_prots.var.index:\n# #     if protein in cite_genes.var.index:\n#     print(protein)\n\n#     y = pd.DataFrame(\n#         cite_prots[:, cite_prots.var.index == protein].X.A, \n#         cite_prots.obs.index, \n#         [protein]\n#     )\n\n# #         df_scores = calc_lgb(X.values, y.values.flatten(), protein,  df_scores)\n#     df_scores = calc_ridge(X.values, y.values.flatten(), protein,  df_scores, 'cell types')","metadata":{"execution":{"iopub.status.busy":"2023-01-11T15:39:32.879858Z","iopub.execute_input":"2023-01-11T15:39:32.880172Z","iopub.status.idle":"2023-01-11T15:39:32.887293Z","shell.execute_reply.started":"2023-01-11T15:39:32.880146Z","shell.execute_reply":"2023-01-11T15:39:32.886386Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Predicting based on - cell type +  own RNA \n\n(own RNA -   RNA corresponding to the protein, e.g. for CD36 we take CD36  )","metadata":{}},{"cell_type":"code","source":"# %%time\n# for protein in cite_prots.var.index:\n#     if protein in cite_genes.var.index:\n#         print(protein)\n#         X = pd.DataFrame(\n#             cite_genes[:, cite_genes.var.index == protein].X.A, \n#             cite_genes.obs.index, \n#             [protein]\n#         )\n#         X = pd.concat([\n#             X, pd.get_dummies(adata.obs[['cell_type']])\n#         ], axis = 1)\n        \n#         y = pd.DataFrame(\n#             cite_prots[:, cite_prots.var.index == protein].X.A, \n#             cite_prots.obs.index, \n#             [protein]\n#         )\n        \n#         df_scores = calc_lgb(X.values, y.values.flatten(), protein,  df_scores, 'CT + ownRNA')\n#         df_scores = calc_ridge(X.values, y.values.flatten(), protein,  df_scores, 'CT + ownRNA')","metadata":{"execution":{"iopub.status.busy":"2023-01-11T15:39:32.888668Z","iopub.execute_input":"2023-01-11T15:39:32.889761Z","iopub.status.idle":"2023-01-11T15:39:32.896554Z","shell.execute_reply.started":"2023-01-11T15:39:32.889700Z","shell.execute_reply":"2023-01-11T15:39:32.894757Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Predicting based on PCA top 100 features","metadata":{}},{"cell_type":"code","source":"# %%time\n\n# sc.pp.pca(cite_genes, n_comps=200) # 7 minutes\n# # 50 comps are already in adata (cite_prots.obsm['GEX_X_pca']), \n# # if you need 50 comps you can commit the above line","metadata":{"execution":{"iopub.status.busy":"2023-01-11T15:39:32.898145Z","iopub.execute_input":"2023-01-11T15:39:32.898748Z","iopub.status.idle":"2023-01-11T15:39:32.907057Z","shell.execute_reply.started":"2023-01-11T15:39:32.898694Z","shell.execute_reply":"2023-01-11T15:39:32.906194Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# X = cite_genes.obsm['X_pca'][:, :100]","metadata":{"execution":{"iopub.status.busy":"2023-01-11T15:39:32.908462Z","iopub.execute_input":"2023-01-11T15:39:32.908991Z","iopub.status.idle":"2023-01-11T15:39:32.915780Z","shell.execute_reply.started":"2023-01-11T15:39:32.908956Z","shell.execute_reply":"2023-01-11T15:39:32.914768Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# %%time\n# for protein in cite_prots.var.index:\n# #     if protein in cite_genes.var.index:\n#     print(protein)\n\n#     y = pd.DataFrame(\n#         cite_prots[:, cite_prots.var.index == protein].X.A, \n#         cite_prots.obs.index, \n#         [protein]\n#     )\n\n#     df_scores = calc_lgb(X, y.values.flatten(), protein,  df_scores, 'pca100')\n#     df_scores = calc_ridge(X, y.values.flatten(), protein,  df_scores, 'pca100')\n","metadata":{"execution":{"iopub.status.busy":"2023-01-11T15:39:32.917351Z","iopub.execute_input":"2023-01-11T15:39:32.917932Z","iopub.status.idle":"2023-01-11T15:39:32.924861Z","shell.execute_reply.started":"2023-01-11T15:39:32.917895Z","shell.execute_reply":"2023-01-11T15:39:32.923977Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"gc.collect()","metadata":{"execution":{"iopub.status.busy":"2023-01-11T15:39:32.926454Z","iopub.execute_input":"2023-01-11T15:39:32.927020Z","iopub.status.idle":"2023-01-11T15:39:33.086749Z","shell.execute_reply.started":"2023-01-11T15:39:32.926987Z","shell.execute_reply":"2023-01-11T15:39:33.085700Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Prediction based on proteins - all except the target one ","metadata":{}},{"cell_type":"code","source":"# %%time\n\n\n# for protein in cite_prots.var.index:\n# #     if protein in cite_genes.var.index:\n#     print(protein)\n\n#     X = pd.DataFrame(\n#         cite_prots.X.A, \n#         cite_prots.obs.index, \n#         cite_prots.var.index.values.flatten()\n#     ).drop(protein , axis = 1)\n\n#     y = pd.DataFrame(\n#         cite_prots[:, cite_prots.var.index == protein].X.A, \n#         cite_prots.obs.index, \n#         [protein]\n#     )\n\n#     df_scores = calc_lgb(X.values, y.values.flatten(), protein,  df_scores, 'allProteins')\n#     df_scores = calc_ridge(X.values, y.values.flatten(), protein,  df_scores, 'allProteins')\n        \n        \n# gc.collect()","metadata":{"execution":{"iopub.status.busy":"2023-01-11T15:39:33.088596Z","iopub.execute_input":"2023-01-11T15:39:33.088930Z","iopub.status.idle":"2023-01-11T15:39:33.095492Z","shell.execute_reply.started":"2023-01-11T15:39:33.088890Z","shell.execute_reply":"2023-01-11T15:39:33.094526Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# l = []\n# for IX in df_scores.index:\n#     if ('Corr' in IX) and ('Spearman' not in IX): l.append(IX)\n# display( df_scores.loc[l,:].sort_values(target_name) )\n","metadata":{"execution":{"iopub.status.busy":"2023-01-11T15:39:33.097253Z","iopub.execute_input":"2023-01-11T15:39:33.097626Z","iopub.status.idle":"2023-01-11T15:39:33.107305Z","shell.execute_reply.started":"2023-01-11T15:39:33.097592Z","shell.execute_reply":"2023-01-11T15:39:33.106286Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Ridge on ALL RNA features ( 11 minutes of calculations ) (catboost)","metadata":{}},{"cell_type":"code","source":"# %%time\n# X = cite_genes.X #.A\n\n# # for protein in cite_prots.var.index:\n# #     if protein in cite_genes.var.index:\n# #         print(protein)\n       \n# #         y = pd.DataFrame(\n# #             cite_prots[:, cite_prots.var.index == protein].X.A, \n# #             cite_prots.obs.index, \n# #             [protein]\n# #         )\n        \n# # #         df_scores = calc_lgb(X, y.values.flatten(), protein,  df_scores, 'allRNA')\n# #         df_scores = calc_ridge(X, y.values.flatten(), protein,  df_scores, 'allRNA')\n# #     gc.collect()\n    \n    \n# !mkdir allrna_models\n\n# cat_imps = pd.DataFrame(index = cite_genes.var.index.tolist())\n# for protein in tqdm(cite_prots.var.index):\n#     if protein in cite_genes.var.index:\n#         print(protein)\n#         cat = CatBoostRegressor(task_type=\"GPU\", devices='0:1')\n       \n#         y = pd.DataFrame(\n#             cite_prots[:, cite_prots.var.index == protein].X.A, \n#             cite_prots.obs.index, \n#             [protein]\n#         )\n        \n# #         df_scores = calc_lgb(X, y.values.flatten(), protein,  df_scores, 'allRNA')\n# #         df_scores = calc_ridge(X, y.values.flatten(), protein,  df_scores, 'allRNA')\n#         cat.fit(X,y,verbose=False, plot=False)\n#         cat.save_model(f'allrna_models/{protein}')\n        \n#         tmp_cat_imps = pd.DataFrame({protein:cat.feature_importances_}, index = cite_genes.var.index.tolist())\n\n#         cat_imps = pd.concat([cat_imps, tmp_cat_imps], axis = 1)\n        \n    \n# #         explainer = shap.Explainer(cat.predict, X)\n# # #         # Calculates the SHAP values - It takes some time\n# # #         shap_values = explainer(X)\n# #         shap_values = explainer.shap_values(X)\n# #         shap.plots.beeswarm(shap_values)\n    \n#     gc.collect()\n#     cat_imps.to_csv('PredictionValuesChange_allRNA.csv')","metadata":{"execution":{"iopub.status.busy":"2023-01-11T15:39:33.108387Z","iopub.execute_input":"2023-01-11T15:39:33.108633Z","iopub.status.idle":"2023-01-11T15:39:33.117828Z","shell.execute_reply.started":"2023-01-11T15:39:33.108606Z","shell.execute_reply":"2023-01-11T15:39:33.116967Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"!mkdir shap_each_cell\nX = cite_genes.X #.A\nshap_imps = pd.DataFrame(index = cite_genes.var.index.tolist())\ncat_imps = pd.DataFrame(index = cite_genes.var.index.tolist())\n\nfor protein in tqdm(cite_prots.var.index):\n    if protein in cite_genes.var.index:\n        print(protein)\n        cat = CatBoostRegressor(task_type=\"GPU\", devices='0:1')\n        cat.load_model(f'/kaggle/input/nips2021-each-prot-predict/allrna_models/{protein}')\n        \n        X = cite_genes.X\n        y = pd.DataFrame(\n            cite_prots[:, cite_prots.var.index == protein].X.A, \n            cite_prots.obs.index, \n            [protein]\n        )\n        \n        tmp_cat_imp = cat.get_feature_importance(Pool(X, y))\n        tmp_cat_imp_df = pd.DataFrame({protein:tmp_cat_imp}, index = cite_genes.var.index.tolist())\n        cat_imps = pd.concat([cat_imps, tmp_cat_imp_df], axis = 1)\n        # del less important genes\n        N = 0.1\n        imp_genes = tmp_cat_imp_df[tmp_cat_imp_df[protein] > N].index.tolist()\n        print(f'number imp genes (PredictionValuesChange > {N}): {len(imp_genes)}')\n        \n        X = cite_genes[:, cite_genes.var.index.isin(imp_genes)].X \n        \n        cat = CatBoostRegressor(task_type=\"GPU\", devices='0:1')\n        cat.fit(X,y,verbose=False, plot=False)\n        \n        # https://github.com/catboost/catboost/blob/master/catboost/tutorials/model_analysis/shap_values_tutorial.ipynb\n        shap_values = cat.get_feature_importance(Pool(X, y), type = 'ShapValues')[:, :-1]\n        \n        pd.DataFrame(shap_values, index = cite_genes.obs.index, columns = imp_genes).to_csv(f'shap_each_cell/{protein}_shaps.csv')\n        tmp_shap_imp_df = pd.DataFrame({protein:shap_values.mean(axis = 0)}, index = imp_genes)\n        shap_imps = pd.concat([shap_imps, tmp_shap_imp_df], axis = 1)\n        \n#         shap.summary_plot(shap_values, X)\n        gc.collect()\nshap_imps.to_csv('shap_vals_allRNA.csv')\ncat_imps.to_csv('PredictionValuesChange_allRNA.csv')","metadata":{"execution":{"iopub.status.busy":"2023-01-11T15:39:33.119522Z","iopub.execute_input":"2023-01-11T15:39:33.119907Z","iopub.status.idle":"2023-01-11T15:41:36.637552Z","shell.execute_reply.started":"2023-01-11T15:39:33.119871Z","shell.execute_reply":"2023-01-11T15:41:36.635668Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# from sklearn.inspection import permutation_importance\n# import shap\n\n# from catboost import CatBoostRegressor, Pool","metadata":{"execution":{"iopub.status.busy":"2023-01-11T15:41:36.639309Z","iopub.status.idle":"2023-01-11T15:41:36.640178Z","shell.execute_reply.started":"2023-01-11T15:41:36.639907Z","shell.execute_reply":"2023-01-11T15:41:36.639934Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# explainer = shap.TreeExplainer(cat)\n# shap_values = explainer.shap_values(X)\n# shap.summary_plot(shap_values, X)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Prediction based on both rna and proteins (except target one ) \n","metadata":{"execution":{"iopub.status.busy":"2022-12-23T18:04:12.858384Z","iopub.execute_input":"2022-12-23T18:04:12.859084Z","iopub.status.idle":"2022-12-23T18:04:12.865728Z","shell.execute_reply.started":"2022-12-23T18:04:12.859047Z","shell.execute_reply":"2022-12-23T18:04:12.86498Z"}}},{"cell_type":"code","source":"# for protein in cite_prots.var.index:\n#     if protein in cite_genes.var.index:\n#         print(protein)\n        \n#         X_prot = pd.DataFrame(\n#             cite_prots.X.A, \n#             cite_prots.obs.index, \n#             cite_prots.var.index.values.flatten()\n#         ).drop(protein , axis = 1)\n        \n#         X_rna = pd.DataFrame(\n#             cite_genes.X.A, \n#             cite_genes.obs.index, \n#             cite_genes.var.index.values.flatten()\n#         ) #.drop(protein , axis = 1)\n        \n#         X = pd.concat( [ X_prot, X_rna], axis = 1)\n#         del X_prot, X_rna \\ gc.collect()\n\n#         y = pd.DataFrame(\n#             cite_prots[:, cite_prots.var.index == protein].X.A, \n#             cite_prots.obs.index, \n#             [protein]\n#         )\n        \n#         df_scores = calc_lgb(X.values, y.values.flatten(), protein,  df_scores, 'allRNA+Proteins')\n#         df_scores = calc_ridge(X.values, y.values.flatten(), protein,  df_scores, 'allRNA+Proteins')\n        \n        \n# gc.collect()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# l = []\n# for IX in df_scores.index:\n#     if ('Corr' in IX) and ('Spearman' not in IX): l.append(IX)\n# display( df_scores.loc[l,:].sort_values(target_name) )\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# df_y.corr()[target_name].sort_values(ascending= False)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Prediction based on RNA (except own RNA) and Proteins (except own proteins) ","metadata":{}},{"cell_type":"code","source":"# %%time\n# for protein in cite_prots.var.index:\n#     if protein in cite_genes.var.index:\n#         print(protein)\n        \n#         X_prot = pd.DataFrame(\n#             cite_prots.X.A, \n#             cite_prots.obs.index, \n#             cite_prots.var.index.values.flatten()\n#         ).drop(protein , axis = 1)\n        \n#         X_rna = pd.DataFrame(\n#             cite_genes.X.A, \n#             cite_genes.obs.index, \n#             cite_genes.var.index.values.flatten()\n#         ).drop(protein , axis = 1)\n        \n#         X = pd.concat( [ X_prot, X_rna], axis = 1)\n#         del X_prot, X_rna \\ gc.collect()\n\n#         y = pd.DataFrame(\n#             cite_prots[:, cite_prots.var.index == protein].X.A, \n#             cite_prots.obs.index, \n#             [protein]\n#         )\n        \n#         df_scores = calc_lgb(X.values, y.values.flatten(), protein,  df_scores, 'Proteins + RNA except own')\n#         df_scores = calc_ridge(X.values, y.values.flatten(), protein,  df_scores, 'Proteins + RNA except own')\n        \n        \n# gc.collect()\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# l = []\n# for IX in df_scores.index:\n#     if ('Corr' in IX) and ('Spearman' not in IX): l.append(IX)\n# display( df_scores.loc[l,:].sort_values(target_name) )","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Prediction based on RNA (all except own RNA) ","metadata":{}},{"cell_type":"code","source":"# import gc\n# gc.collect()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# %%time\n\n# for protein in cite_prots.var.index:\n#     if protein in cite_genes.var.index:\n#         print(protein)\n#         X = pd.DataFrame(\n#             cite_genes.X.A, \n#             cite_genes.obs.index, \n#             cite_genes.var.index.values.flatten()\n#         ).drop(protein , axis = 1)\n       \n#         y = pd.DataFrame(\n#             cite_prots[:, cite_prots.var.index == protein].X.A, \n#             cite_prots.obs.index, \n#             [protein]\n#         )\n        \n# #         df_scores = calc_lgb(X, y.values.flatten(), protein,  df_scores, 'all RNA except own')\n#         df_scores = calc_ridge(X, y.values.flatten(), protein,  df_scores, 'all RNA except own')\n#     gc.collect()","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# l = []\n# for IX in df_scores.index:\n#     if ('Corr' in IX) and ('Spearman' not in IX): l.append(IX)\n# display( df_scores.loc[l,:].sort_values(target_name) )","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Save resuls","metadata":{"execution":{"iopub.status.busy":"2022-12-23T17:06:55.028528Z","iopub.execute_input":"2022-12-23T17:06:55.02897Z","iopub.status.idle":"2022-12-23T17:06:55.034005Z","shell.execute_reply.started":"2022-12-23T17:06:55.028933Z","shell.execute_reply":"2022-12-23T17:06:55.032823Z"}}},{"cell_type":"code","source":"# df_scores.to_csv('Scores_all_prots.csv')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# df = pd.read_csv('/kaggle/input/simple-models-for-cds-on-nips21/Scores_all_prots.csv', index_col = 0).iloc[2:, :]","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# for protein in df.columns:\n# #     break\n#     print(f'Top 5 scores for {protein}: \\n {df[protein].sort_values(ascending = False).iloc[0:5]} \\n')","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"t = time.time()-t0start\nprint('%.1f hours  = %.1f minutes = %.1f seconds passed total '%( t/3600,t/60, t ) )","metadata":{"trusted":true},"execution_count":null,"outputs":[]}]}