{"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\nDimensional reduction and visualization proteins features for CAFA5 Kaggle challenge - https://www.kaggle.com/competitions/cafa-5-protein-function-prediction\n    \nThe present notebook uses protein embeddings calculated by esm2 protein language model the Meta. \nComputed and kindly shared with Kaggle community by Andrey: https://www.kaggle.com/competitions/cafa-5-protein-function-prediction/discussion/406168\nplease upvote his post and datasets ! \n\n\n\nPS\n\nSee also similar script for the other emebddings: \n\nSGT https://www.kaggle.com/alexandervc/cafa5-sgt-embeddings-visualizations\n\nLevenstein distance features: https://www.kaggle.com/code/alexandervc/cafa5-levenshtein-distances-feature-visualizations\n\nsame esm2 with TSNE by Tilli: https://www.kaggle.com/competitions/cafa-5-protein-function-prediction/discussion/406336\n\nT5 https://www.kaggle.com/alexandervc/cafa5-t5-protein-embeddings-visualization\n\nProtBert: https://www.kaggle.com/alexandervc/cafa5-protbert-protein-embeddings-visualization\n\nT5, ProtBert, etc: https://www.kaggle.com/code/alexandervc/cafa5-visualization-dimensional-reductions-embed\n\n\n\nPS\n\nHuge collection of various dimensional reduction / visualization methods : \n\nhttps://www.kaggle.com/code/alexandervc/cite-seq2023-dim-red-visualizations\n\nhttps://www.kaggle.com/code/alexandervc/medulloblastoma-gse85217cavalli-dimred-visualiz\n","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\n\nimport time\nt0start = time.time() \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 'png' not in 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","execution":{"iopub.status.busy":"2023-05-06T10:04:51.50358Z","iopub.execute_input":"2023-05-06T10:04:51.504057Z","iopub.status.idle":"2023-05-06T10:04:51.653133Z","shell.execute_reply.started":"2023-05-06T10:04:51.504013Z","shell.execute_reply":"2023-05-06T10:04:51.651901Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import matplotlib.pyplot as plt\nimport seaborn as sns","metadata":{"execution":{"iopub.status.busy":"2023-05-02T08:24:04.083183Z","iopub.execute_input":"2023-05-02T08:24:04.083585Z","iopub.status.idle":"2023-05-02T08:24:04.08995Z","shell.execute_reply.started":"2023-05-02T08:24:04.08355Z","shell.execute_reply":"2023-05-02T08:24:04.088591Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Load pre-computed embeddings data","metadata":{}},{"cell_type":"code","source":"%%time\n\nfeatures = 'ProtBert' # 'T5' #  'esm2_t33_650M'\n\nif features == 'T5':  \n    fn = '/kaggle/input/t5embeds/train_embeds.npy'\n    fn4submit = '/kaggle/input/t5embeds/test_embeds.npy'\n    str_data_inf = 'CAFA5 T5 embeddings '\n    \nelif features == 'ProtBert':\n    \n    fn = '/kaggle/input/protbert-embeddings-for-cafa5/train_embeddings.npy'\n    fn4submit = '/kaggle/input/protbert-embeddings-for-cafa5/test_embeddings.npy'\n    str_data_inf = 'CAFA5 ProtBert embeddings '\n\nelif features == 'esm2_t33_650M':\n    fn =        '/kaggle/input/23468234/train_embeds_esm2_t33_650M_UR50D.npy'\n    fn4submit = '/kaggle/input/23468234/test_embeds_esm2_t33_650M_UR50D.npy'\n    str_data_inf = 'CAFA5 ' + features\n\n    # /kaggle/input/23468234/train_ids_esm2_t33_650M_UR50D.npy\n    # /kaggle/input/23468234/test_ids_esm2_t33_650M_UR50D.npy\n    # /kaggle/input/23468234/test_embeds_esm2_t33_650M_UR50D.npy\n    # /kaggle/input/23468234/train_embeds_esm2_t33_650M_UR50D.npy\n\nelif features == 'esm2_t30_150M':\n    fn = '/kaggle/input/8947923/train_embeds_esm2_t30_150M_UR50D.npy'\n    fn4submit = '/kaggle/input/8947923/test_embeds_esm2_t30_150M_UR50D.npy'\n    str_data_inf = 'CAFA5 ' + features\n\n    # /kaggle/input/8947923/train_embeds_esm2_t30_150M_UR50D.npy\n    # /kaggle/input/8947923/test_embeds_esm2_t30_150M_UR50D.npy\n\n    # /kaggle/input/8947923/train_ids_esm2_t30_150M_UR50D.npy\n    # /kaggle/input/8947923/test_ids_esm2_t30_150M_UR50D.npy\n\nelif features == 'esm2_t12_35M':\n    fn =        '/kaggle/input/3023750/train_embeds_esm2_t12_35M_UR50D.npy'\n    fn4submit = '/kaggle/input/3023750/test_embeds_esm2_t12_35M_UR50D.npy'\n    str_data_inf = 'CAFA5 ' + features\n    # /kaggle/input/3023750/train_embeds_esm2_t12_35M_UR50D.npy\n    # /kaggle/input/3023750/test_embeds_esm2_t12_35M_UR50D.npy\n    # /kaggle/input/3023750/test_ids_esm2_t12_35M_UR50D.npy\n    # /kaggle/input/3023750/train_ids_esm2_t12_35M_UR50D.npy\n\nelif features == 'esm2_t6_8M':\n    fn = '/kaggle/input/315701375/train_embeds_esm2_t6_8M_UR50D.npy'\n    fn4submit = '/kaggle/input/315701375/test_embeds_esm2_t6_8M_UR50D.npy'\n    str_data_inf = 'CAFA5 ' + features\n    # /kaggle/input/315701375/train_embeds_esm2_t6_8M_UR50D.npy\n    # /kaggle/input/315701375/test_embeds_esm2_t6_8M_UR50D.npy\n    # /kaggle/input/315701375/train_ids_esm2_t6_8M_UR50D.npy        \n    # /kaggle/input/315701375/test_ids_esm2_t6_8M_UR50D.npy\n\nelif features == 'SGTv12':\n    fn = '/kaggle/input/cafa5-sgt-protein-embeddings/run12_fit_on_15000x2/sgt_embendings_train_15000_15000.csv' \n    fn4submit = '/kaggle/input/cafa5-sgt-protein-embeddings/run12_fit_on_15000x2/sgt_embendings_test_15000_15000.csv'\n    str_data_inf = 'CAFA5 SGT embeddings v12 '\n\nelif features == 'Leven5000':\n    fn = '/kaggle/input/cafa5-levenshtein-distance-features-big/run9_achnor_set5000/df_Levenshtein_distance_5000_features_train.csv' \n    fn4submit = '/kaggle/input/cafa5-levenshtein-distance-features-big/run9_achnor_set5000/df_Levenshtein_distance_5000_features_test.csv'\n    str_data_inf = 'CAFA5 Levenshtein distance features anchor 5000 '\n\nelif features == 'Leven3000':\n    fn = '/kaggle/input/cafa5-levenshtein-distance-features-big/run8_anchor_set_3000/df_Levenshtein_distance_3000_features_train.csv' \n    fn4submit = '/kaggle/input/cafa5-levenshtein-distance-features-big/run8_anchor_set_3000/df_Levenshtein_distance_3000_features_test.csv'\n    str_data_inf = 'CAFA5 Levenshtein distance features anchor 3000 '\n\nelif features == 'Leven1000':\n    fn = '/kaggle/input/cafa5-features-etc/df_Levenshtein_distance_1000_features_train.csv' \n    fn4submit = '/kaggle/input/cafa5-features-etc/df_Levenshtein_distance_1000_features_test.csv'\n    str_data_inf = 'CAFA5 Levenshtein distance features anchor 1000 '\n\nelif features == 'Leven100':\n    fn = '/kaggle/input/cafa5-features-etc/df_Levenshtein_distance_100_features_train.csv' \n    fn4submit = '/kaggle/input/cafa5-features-etc/df_Levenshtein_distance_100_features_test.csv'\n    str_data_inf = 'CAFA5 Levenshtein distance features anchor 100 '\n# else:\n#     fn = '/kaggle/input/t5embeds/train_embeds.npy'\n#     features = 'T5'\n    \n\ndef load_features(fn ):\n    print(fn)\n    if '.csv' in fn:\n        df = pd.read_csv(fn, index_col = 0)\n    #     df = df.rank(axis = 1 )\n        X = df.values\n    elif '.npy' in fn:\n        X = np.load(fn)\n    print(X.shape)\n    return X\n\nX = load_features(fn)\nstr_data_inf += ' train'\nprint(X.shape)\nX[:10,:10]","metadata":{"execution":{"iopub.status.busy":"2023-05-02T08:24:19.780061Z","iopub.execute_input":"2023-05-02T08:24:19.780466Z","iopub.status.idle":"2023-05-02T08:24:20.315993Z","shell.execute_reply.started":"2023-05-02T08:24:19.780432Z","shell.execute_reply":"2023-05-02T08:24:20.315001Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# PCA","metadata":{}},{"cell_type":"code","source":"import matplotlib.pyplot as plt\nimport seaborn as sns\n\nmarker1 =  'o' # 'o' #\npalette1 = 'rainbow' #  \"viridis\"  #  'rainbow' #  \"viridis\"# \nalpha1 = 1 ","metadata":{"execution":{"iopub.status.busy":"2023-05-01T15:07:20.639465Z","iopub.execute_input":"2023-05-01T15:07:20.640078Z","iopub.status.idle":"2023-05-01T15:07:20.646425Z","shell.execute_reply.started":"2023-05-01T15:07:20.640027Z","shell.execute_reply":"2023-05-01T15:07:20.645001Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\nfrom sklearn.decomposition import PCA\npca = PCA(n_components=6)\nr = pca.fit_transform(X)\n\nprint(np.sum(pca.explained_variance_ratio_))\n","metadata":{"execution":{"iopub.status.busy":"2023-05-01T15:07:20.6498Z","iopub.execute_input":"2023-05-01T15:07:20.650283Z","iopub.status.idle":"2023-05-01T15:07:40.841735Z","shell.execute_reply.started":"2023-05-01T15:07:20.650233Z","shell.execute_reply":"2023-05-01T15:07:40.840399Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"v4color = pd.Series(r[:,0])\nv4color.name = 'PCA1'","metadata":{"execution":{"iopub.status.busy":"2023-05-01T15:07:40.84353Z","iopub.execute_input":"2023-05-01T15:07:40.844282Z","iopub.status.idle":"2023-05-01T15:07:40.850849Z","shell.execute_reply.started":"2023-05-01T15:07:40.844234Z","shell.execute_reply":"2023-05-01T15:07:40.849873Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time \nn_x_subplots = 2\nc = 0\ncc = 0\nfor (i,j) in [(0,1),(0,2),(0,3),(1,2),(1,3),(2,3),(0,4),(1,4),(2,4),(3,4), (0,5)]:\n    cc += 1\n    if c % n_x_subplots == 0:\n        if c > 0:\n            plt.show()\n        fig = plt.figure(figsize = (20,5) ); c = 0\n        plt.suptitle(str_data_inf +  ' n_samples='+str(len(r)),fontsize = 20 ) \n\n\n    c += 1; fig.add_subplot(1,n_x_subplots ,c)\n\n    ax = sns.scatterplot(x = r[:,i],y=r[:,j], hue = v4color,  palette = palette1,  marker = marker1 , alpha = alpha1 )\n    if c == 1:\n        plt.setp(ax.get_legend().get_texts(), fontsize=12) # for legend text\n        plt.setp(ax.get_legend().get_title(), fontsize=12) # for legend title    \n        plt.legend( bbox_to_anchor = (1.15, 1.), loc='upper right')\n    else:\n        plt.legend('')\n    plt.title('PCA '+str(i)+' '+str(j), fontsize = 20)\n    \n    \n    \nplt.show() ","metadata":{"execution":{"iopub.status.busy":"2023-05-01T15:07:40.853078Z","iopub.execute_input":"2023-05-01T15:07:40.853957Z","iopub.status.idle":"2023-05-01T15:08:57.605954Z","shell.execute_reply.started":"2023-05-01T15:07:40.853908Z","shell.execute_reply":"2023-05-01T15:08:57.60494Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# ICA","metadata":{}},{"cell_type":"code","source":"%%time\nfrom sklearn.decomposition import FastICA\nreducer = FastICA(n_components=5, random_state=0, whiten='unit-variance')\n\nr = reducer.fit_transform(X)\n\nn_x_subplots = 2\nc = 0\ncc = 0\n\nfor (i,j) in [(0,1),(0,2),(1,2),(2,3),(2,4),(3,4)]:\n\n    if c % n_x_subplots == 0:\n        if c > 0:\n            plt.show()\n        fig = plt.figure(figsize = (20,5) ); c = 0\n        plt.suptitle(str_data_inf + ' n_samples='+str(len(r)),fontsize = 20 ) \n\n\n    c += 1; fig.add_subplot(1,n_x_subplots ,c)\n\n    ax = sns.scatterplot(x = r[:,i],y=r[:,j], hue = v4color,  palette = palette1, marker =  marker1, alpha = alpha1  )\n    if 1:\n        plt.setp(ax.get_legend().get_texts(), fontsize=12) # for legend text\n        plt.setp(ax.get_legend().get_title(), fontsize=12) # for legend title        \n    plt.title('ICA '+str(i)+' '+str(j), fontsize = 20)\n    if c == 1:\n        plt.setp(ax.get_legend().get_texts(), fontsize=12) # for legend text\n        plt.setp(ax.get_legend().get_title(), fontsize=12) # for legend title    \n        plt.legend( bbox_to_anchor = (1.15, 1.), loc='upper right')\n    else:\n        plt.legend('')\n    \nplt.show() ","metadata":{"execution":{"iopub.status.busy":"2023-05-01T09:07:33.243305Z","iopub.execute_input":"2023-05-01T09:07:33.243666Z","iopub.status.idle":"2023-05-01T09:14:42.987828Z","shell.execute_reply.started":"2023-05-01T09:07:33.243631Z","shell.execute_reply":"2023-05-01T09:14:42.986682Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# UMAP ","metadata":{}},{"cell_type":"code","source":"%%time\nimport matplotlib.pyplot as plt\nimport seaborn as sns\n\nimport umap \n\nr = umap.UMAP().fit_transform(X)\n\n# v4color = df_corpus['sequence'].apply(lambda x: len(x))\n# v4color.name = 'sequence length'\n\nprint(r.shape)\n\nfor i,j in [(0,1)]: # ,(0,2),(0,3),(1,2),(1,3),(2,3)]:\n    sns.scatterplot(x = r[:,i], y = r[:,j], hue = v4color, palette = 'rainbow' )\n    plt.xlabel('UMAP'+str(i+1), fontsize = 20)\n    plt.ylabel('UMAP'+str(j+1), fontsize = 20)\n    plt.show()","metadata":{"execution":{"iopub.status.busy":"2023-05-01T09:14:42.989433Z","iopub.execute_input":"2023-05-01T09:14:42.990486Z","iopub.status.idle":"2023-05-01T09:22:46.62777Z","shell.execute_reply.started":"2023-05-01T09:14:42.990445Z","shell.execute_reply":"2023-05-01T09:22:46.626253Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# NCVis - similar to TSNE, bust much faster (faster than UMAP also). From Skolkovo team\n\nhttps://arxiv.org/abs/2001.11411\n\nAleksandr Artemenkov, Maxim Panov\n\nModern methods for data visualization via dimensionality reduction, such as t-SNE, usually have performance issues that prohibit their application to large amounts of high-dimensional data. In this work, we propose NCVis -- a high-performance dimensionality reduction method built on a sound statistical basis of noise contrastive estimation. We show that NCVis outperforms state-of-the-art techniques in terms of speed while preserving the representation quality of other methods. In particular, the proposed approach successfully proceeds a large dataset of more than 1 million news headlines in several minutes and presents the underlying structure in a human-readable way. Moreover, it provides results consistent with classical methods like t-SNE on more straightforward datasets like images of hand-written digits. We believe that the broader usage of such software can significantly simplify the large-scale data analysis and lower the entry barrier to this area.","metadata":{}},{"cell_type":"code","source":"%%time\n!pip install ncvis\nimport ncvis\n\nreducer = ncvis.NCVis()\nr = reducer.fit_transform(X)\n\nfor (i,j) in [(0,1)]:#,(0,2),(1,2),(2,3),(2,4),(3,4)]:\n    fig = plt.figure(figsize = (20,12) ); c = 0\n    ax = sns.scatterplot(x = r[:,i],y=r[:,j], hue = v4color,  palette = palette1, marker =  marker1 , alpha = alpha1  )\n    if 1:\n        plt.setp(ax.get_legend().get_texts(), fontsize=20) # for legend text\n        plt.setp(ax.get_legend().get_title(), fontsize=20) # for legend title        \n    plt.title(' ncvis ' + str_data_inf + ' n_samples='+str(len(r)) , fontsize = 20)\n    plt.show()","metadata":{"execution":{"iopub.status.busy":"2023-05-01T13:45:05.381411Z","iopub.execute_input":"2023-05-01T13:45:05.381949Z","iopub.status.idle":"2023-05-01T13:48:52.141561Z","shell.execute_reply.started":"2023-05-01T13:45:05.381901Z","shell.execute_reply":"2023-05-01T13:48:52.140193Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# pacmap\n\nhttps://arxiv.org/abs/2012.04456\n\nUnderstanding How Dimension Reduction Tools Work: An Empirical Approach to Deciphering t-SNE, UMAP, TriMAP, and PaCMAP for Data Visualization\nYingfan Wang, Haiyang Huang, Cynthia Rudin, Yaron Shaposhnik\nDimension reduction (DR) techniques such as t-SNE, UMAP, and TriMAP have demonstrated impressive visualization performance on many real world datasets. One tension that has always faced these methods is the trade-off between preservation of global structure and preservation of local structure: these methods can either handle one or the other, but not both. In this work, our main goal is to understand what aspects of DR methods are important for preserving both local and global structure: it is difficult to design a better method without a true understanding of the choices we make in our algorithms and their empirical impact on the lower-dimensional embeddings they produce. Towards the goal of local structure preservation, we provide several useful design principles for DR loss functions based on our new understanding of the mechanisms behind successful DR methods. Towards the goal of global structure preservation, our analysis illuminates that the choice of which components to preserve is important. We leverage these insights to design a new algorithm for DR, called Pairwise Controlled Manifold Approximation Projection (PaCMAP), which preserves both local and global structure. Our work provides several unexpected insights into what design choices both to make and avoid when constructing DR algorithms.","metadata":{}},{"cell_type":"code","source":"%%time\n!pip install pacmap\n\nimport pacmap\n\nreducer = pacmap.PaCMAP()\nr = reducer.fit_transform(X)\n\nfor (i,j) in [(0,1)]:#,(0,2),(1,2),(2,3),(2,4),(3,4)]:\n    fig = plt.figure(figsize = (20,12) ); c = 0\n    ax = sns.scatterplot(x = r[:,i],y=r[:,j], hue = v4color,  palette = palette1, marker =  marker1 , alpha = alpha1 )\n    if 1:\n        plt.setp(ax.get_legend().get_texts(), fontsize=20) # for legend text\n        plt.setp(ax.get_legend().get_title(), fontsize=20) # for legend title        \n    plt.title(' pacmap ' + str_data_inf  + ' n_samples='+str(len(r)), fontsize = 20)\n    plt.show()","metadata":{"execution":{"iopub.status.busy":"2023-05-01T09:26:24.277717Z","iopub.execute_input":"2023-05-01T09:26:24.278814Z","iopub.status.idle":"2023-05-01T09:30:06.04836Z","shell.execute_reply.started":"2023-05-01T09:26:24.278764Z","shell.execute_reply":"2023-05-01T09:30:06.04676Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Trimap","metadata":{}},{"cell_type":"code","source":"%%time\n!pip install trimap \nimport trimap\n\nreducer = trimap.TRIMAP()\nr = reducer.fit_transform(X)\n\nfor (i,j) in [(0,1)]:#,(0,2),(1,2),(2,3),(2,4),(3,4)]:\n    fig = plt.figure(figsize = (20,12) ); c = 0\n    ax = sns.scatterplot(x = r[:,i],y=r[:,j], hue = v4color,  palette = palette1, marker =  marker1 , alpha = alpha1 )\n    if 1:\n        plt.setp(ax.get_legend().get_texts(), fontsize=20) # for legend text\n        plt.setp(ax.get_legend().get_title(), fontsize=20) # for legend title        \n    plt.title(' trimap ' + str_data_inf  + ' n_samples='+str(len(r)), fontsize = 20)\n    plt.show()","metadata":{"execution":{"iopub.status.busy":"2023-05-01T09:30:06.050194Z","iopub.execute_input":"2023-05-01T09:30:06.050598Z","iopub.status.idle":"2023-05-01T09:35:03.986541Z","shell.execute_reply.started":"2023-05-01T09:30:06.05056Z","shell.execute_reply":"2023-05-01T09:35:03.98487Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# MulticoreTSNE","metadata":{}},{"cell_type":"code","source":"%%time\n!pip install MulticoreTSNE\nfrom MulticoreTSNE import MulticoreTSNE as TSNE\n\nreducer = TSNE(n_jobs=4)\nr = reducer.fit_transform(X)\n\nfor (i,j) in [(0,1)]:#,(0,2),(1,2),(2,3),(2,4),(3,4)]:\n    fig = plt.figure(figsize = (20,12) ); c = 0\n    ax = sns.scatterplot(x = r[:,i],y=r[:,j], hue = v4color,  palette = palette1 , marker =  marker1 , alpha = alpha1 )\n    if 1:\n        plt.setp(ax.get_legend().get_texts(), fontsize=20) # for legend text\n        plt.setp(ax.get_legend().get_title(), fontsize=20) # for legend title        \n    \n    plt.title(' MulticoreTSNE ' + str_data_inf +  ' n_samples='+str(len(r)), fontsize = 20)\n    plt.show()","metadata":{"execution":{"iopub.status.busy":"2023-05-01T09:35:03.99264Z","iopub.execute_input":"2023-05-01T09:35:03.993224Z","iopub.status.idle":"2023-05-01T09:35:20.50432Z","shell.execute_reply.started":"2023-05-01T09:35:03.993165Z","shell.execute_reply":"2023-05-01T09:35:20.502927Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# UMAP several n_neighbors min_dis","metadata":{}},{"cell_type":"code","source":"%%time\n# 8 min for (70988, 140) - six times \nfig = plt.figure(figsize = (20,16)); c = 0; cc = 0\nplt.suptitle('UMAP ' + str_data_inf+   ' n_samples='+str(len(X)) , fontsize = 20 )\nfor min_dist in [0.1, 0.9]: # 0.1 - defailt min_dist\n    for n_neighbors in [5,15, 100]: # 15 - default n_neighbors,\n        c += 1; fig.add_subplot(2,3,c)\n        str_inf = 'n_neighbors='+str(n_neighbors) + ' min_dist='+str(min_dist) \n        reducer = umap.UMAP(n_neighbors = n_neighbors, min_dist = min_dist, n_components= 2  )# random_state=42 # metric = \n        r2 = reducer.fit_transform(X)\n\n        i,j = 0,1\n        ax = sns.scatterplot(x = r2[:,i],y = r2[:,j], hue = v4color, palette = palette1 , marker =  marker1 , alpha = alpha1 )\n        if 1:\n            plt.setp(ax.get_legend().get_texts(), fontsize=12) # for legend text\n            plt.setp(ax.get_legend().get_title(), fontsize=12) # for legend title        \n        plt.title(str_inf , fontsize = 20)\n        plt.xlabel('UMAP'+str(i+1), fontsize = 20)\n        plt.ylabel('UMAP'+str(j+1), fontsize = 20)\n        \n#         if cc>0:\n#             plt.legend('')\n#         cc+=1\n        if c == 1:\n            plt.setp(ax.get_legend().get_texts(), fontsize=10) # for legend text\n            plt.setp(ax.get_legend().get_title(), fontsize=10) # for legend title    \n            plt.legend( bbox_to_anchor = (1.15, 1.), loc='upper right')\n        else:\n            plt.legend('')\n\n\nplt.show() ","metadata":{"execution":{"iopub.status.busy":"2023-05-01T09:35:20.505928Z","iopub.execute_input":"2023-05-01T09:35:20.506307Z","iopub.status.idle":"2023-05-01T11:21:17.950224Z","shell.execute_reply.started":"2023-05-01T09:35:20.506273Z","shell.execute_reply":"2023-05-01T11:21:17.94696Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# UMAP dense map ","metadata":{}},{"cell_type":"code","source":"%%time\nimport umap\nreducer = umap.UMAP(densmap=True, random_state=42)\nr2 = reducer.fit_transform(X)\nprint(X.shape)\n\nfor (i,j) in [(0,1)]:# ,(0,2),(1,2), (0,3),(1,3),(2,3), (0,4),(1,4),(2,4),(3,4)]:\n    plt.figure(figsize = (20,10))\n    ax = sns.scatterplot(x = r2[:,i],y = r2[:,j], hue = v4color, palette = palette1 , marker =  marker1 , alpha = alpha1 )\n    if 1:\n        plt.setp(ax.get_legend().get_texts(), fontsize=20) # for legend text\n        plt.setp(ax.get_legend().get_title(), fontsize=20) # for legend title        \n    plt.title('UMAP densmap On ' + str_data_inf +   ' n_samples='+str(len(X)), fontsize = 20)\n    plt.xlabel('UMAP'+str(i+1), fontsize = 12)\n    plt.ylabel('UMAP'+str(j+1), fontsize = 12)\n    plt.show() ","metadata":{"execution":{"iopub.status.busy":"2023-05-01T11:21:17.95342Z","iopub.execute_input":"2023-05-01T11:21:17.953976Z","iopub.status.idle":"2023-05-01T11:48:13.163513Z","shell.execute_reply.started":"2023-05-01T11:21:17.953928Z","shell.execute_reply":"2023-05-01T11:48:13.159921Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"#  UMAP densemap several n_neighbors min_dis","metadata":{}},{"cell_type":"code","source":"%%time\n# 8 min for (70988, 140) - six times \nfig = plt.figure(figsize = (20,16)); c = 0; cc = 0\nplt.suptitle('UMAP densmap ' + str_data_inf+   ' n_samples='+str(len(X)) , fontsize = 20 )\nfor min_dist in [0.1, 0.9]: # 0.1 - defailt min_dist\n    for n_neighbors in [5,15, 100]: # 15 - default n_neighbors,\n        c += 1; fig.add_subplot(2,3,c)\n        str_inf = 'n_neighbors='+str(n_neighbors) + ' min_dist='+str(min_dist) \n        reducer = umap.UMAP(densmap=True, n_neighbors = n_neighbors, min_dist = min_dist, n_components= 2  )# random_state=42 # metric = \n        r2 = reducer.fit_transform(X)\n\n        i,j = 0,1\n        ax = sns.scatterplot(x = r2[:,i],y = r2[:,j], hue = v4color, palette = palette1 , marker =  marker1 , alpha = alpha1 )\n        if 1:\n            plt.setp(ax.get_legend().get_texts(), fontsize=12) # for legend text\n            plt.setp(ax.get_legend().get_title(), fontsize=12) # for legend title        \n        plt.title(str_inf , fontsize = 20)\n        plt.xlabel('UMAP'+str(i+1), fontsize = 20)\n        plt.ylabel('UMAP'+str(j+1), fontsize = 20)\n        \n#         if cc>0:\n#             plt.legend('')\n#         cc+=1\n        if c == 1:\n            plt.setp(ax.get_legend().get_texts(), fontsize=10) # for legend text\n            plt.setp(ax.get_legend().get_title(), fontsize=10) # for legend title    \n            plt.legend( bbox_to_anchor = (1.15, 1.), loc='upper right')\n        else:\n            plt.legend('')\n\n\nplt.show() ","metadata":{"execution":{"iopub.status.busy":"2023-05-01T11:48:13.168536Z","iopub.execute_input":"2023-05-01T11:48:13.168955Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# OpenTSNE","metadata":{}},{"cell_type":"code","source":"%%time\n!pip install opentsne\nfrom openTSNE import TSNE\n\nif str_data_inf != 'Medulloblastoma Intergrated GSE124814': # Fail in that case by unclear reason \n    from openTSNE import TSNE\n\n    reducer = TSNE(\n        perplexity=30,\n        metric=\"euclidean\",\n        n_jobs=8,\n        random_state=42,\n        verbose=True,\n    )\n    r = reducer.fit(X)\n    for (i,j) in [(0,1)]:#,(0,2),(1,2),(2,3),(2,4),(3,4)]:\n        fig = plt.figure(figsize = (20,12) ); c = 0\n        ax = sns.scatterplot(x = r[:,i],y=r[:,j], hue = v4color,  palette = palette1 , marker =  marker1 , alpha = alpha1 )\n        if 1:\n            plt.setp(ax.get_legend().get_texts(), fontsize=12) # for legend text\n            plt.setp(ax.get_legend().get_title(), fontsize=12) # for legend title        \n        plt.title(' openTSNE '+ str_data_inf + ' n_samples='+str(len(r)), fontsize = 20)\n        plt.show()","metadata":{"execution":{"iopub.status.busy":"2023-05-01T13:48:52.144575Z","iopub.execute_input":"2023-05-01T13:48:52.145055Z","iopub.status.idle":"2023-05-01T14:40:18.510019Z","shell.execute_reply.started":"2023-05-01T13:48:52.14501Z","shell.execute_reply":"2023-05-01T14:40:18.508641Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"#  Several other fast methods from sklearn ","metadata":{}},{"cell_type":"code","source":"%%time\n# Based on: \n# https://scikit-learn.org/stable/auto_examples/manifold/plot_compare_methods.html#sphx-glr-auto-examples-manifold-plot-compare-methods-py\n# See also:\n# https://scikit-learn.org/stable/auto_examples/manifold/plot_lle_digits.html#\n\n\n\n# To speed-up reduce dimensions by PCA first\n# X_save = X.copy( )\n#r = pca.fit_transform(X)\n#X = r[:1000,:20]\n\n\nimport time\n\nimport umap \nfrom sklearn import manifold\nfrom sklearn.decomposition import PCA\nfrom sklearn.decomposition import FactorAnalysis\nfrom sklearn.decomposition import NMF\nfrom sklearn.decomposition import FastICA\nfrom sklearn.decomposition import FactorAnalysis\nfrom sklearn.decomposition import LatentDirichletAllocation\nfrom sklearn.ensemble import RandomTreesEmbedding\nfrom sklearn.random_projection import SparseRandomProjection\nfrom sklearn.discriminant_analysis import LinearDiscriminantAnalysis\n\nfrom sklearn.pipeline import make_pipeline\nfrom sklearn.decomposition import TruncatedSVD\n\n\nfrom collections import OrderedDict\nfrom functools import partial\nfrom matplotlib.ticker import NullFormatter\n\n\nn_neighbors = 10\nn_components = 2\n# Set-up manifold methods\nLLE = partial(manifold.LocallyLinearEmbedding,\n              n_neighbors, n_components, eigen_solver='auto')\n\nmethods = OrderedDict()\nmethods['PCA'] = PCA()\nmethods['umap'] = umap.UMAP(n_components = n_components)\nmethods['t-SNE'] = manifold.TSNE(n_components=n_components, init='pca', random_state=0)\nmethods['ICA'] = FastICA(n_components=n_components,         random_state=0)\nmethods['FA'] = FactorAnalysis(n_components=n_components, random_state=0)\n#methods['LLE'] = LLE(method='standard')\n#methods['Modified LLE'] = LLE(method='modified')\n#methods['Isomap'] = manifold.Isomap(n_neighbors, n_components)\nmethods['MDS'] = manifold.MDS(n_components, max_iter=100, n_init=1)\nmethods['SE'] = manifold.SpectralEmbedding(n_components=n_components,\n                                           n_neighbors=n_neighbors)\nmethods['NMF'] = NMF(n_components=n_components,  init='random', random_state=0) \nmethods['RandProj'] = SparseRandomProjection(n_components=n_components, random_state=42)\n\nrand_trees_embed = make_pipeline(RandomTreesEmbedding(n_estimators=200, random_state=0, max_depth=5), TruncatedSVD(n_components=n_components) )\nmethods['RandTrees'] = rand_trees_embed\nmethods['LatDirAll'] = LatentDirichletAllocation(n_components=n_components,  random_state=0)\n#methods['LTSA'] = LLE(method='ltsa') \n#methods['Hessian LLE'] = LLE(method='hessian') \n\nlist_fast_methods = ['FA', 'RandProj','RandTrees',] #  ['PCA','umap','FA', 'NMF','RandProj','RandTrees','ICA'] # 'ICA',\nlist_slow_methods = ['t-SNE','LLE','Modified LLE','Isomap','MDS','SE','LatDirAll','LTSA','Hessian LLE']\n\n# transformer = NeighborhoodComponentsAnalysis(init='random',  n_components=2, random_state=0) # Cannot be applied since supervised - requires y \n# methods['LinDisA'] = LinearDiscriminantAnalysis(n_components=n_components)# Cannot be applied since supervised - requires y \n\n\n# Create figure\nfig = plt.figure(figsize=(25, 16))\nplt.suptitle(str_data_inf + ' n_samples='+str(len(X)) , fontsize = 20 ) \n\n# Plot results\nc = 0\nfor i, (label, method) in enumerate(methods.items()):\n    if label not in  list_fast_methods: #  list_slow_methods :\n        continue\n        \n    t0 = time.time()\n    try:\n        r = method.fit_transform(X)\n    except:\n        print('Got Exception', label )\n        continue \n    t1 = time.time()\n    print(\"%s: %.2g sec\" % (label, t1 - t0))\n    c+=1\n    fig.add_subplot(2, 4 , c) \n    sns.scatterplot(x=r[:,0], y=r[:,1] , hue =  v4color,  palette = palette1, marker =  marker1 , alpha = alpha1 )\n    plt.title(label ,fontsize = 20 )\n    #plt.legend('')\n    if c == 1:\n        plt.setp(ax.get_legend().get_texts(), fontsize=12) # for legend text\n        plt.setp(ax.get_legend().get_title(), fontsize=12) # for legend title    \n        plt.legend( bbox_to_anchor = (1.15, 1.), loc='upper right')\n    else:\n        plt.legend('')\n    \n\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-05-01T14:40:18.51198Z","iopub.execute_input":"2023-05-01T14:40:18.512354Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"#  Several other slow methods from sklearn ","metadata":{}},{"cell_type":"code","source":"N = 2000 # consider only first N samples otherwise methods are tooooo slow","metadata":{"execution":{"iopub.status.busy":"2023-05-01T15:08:57.607385Z","iopub.execute_input":"2023-05-01T15:08:57.607985Z","iopub.status.idle":"2023-05-01T15:08:57.61314Z","shell.execute_reply.started":"2023-05-01T15:08:57.607949Z","shell.execute_reply":"2023-05-01T15:08:57.611741Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## TSNE , LatDirAlloc , Spectral Embedding","metadata":{}},{"cell_type":"code","source":"%%time\n\n# Based on: \n# https://scikit-learn.org/stable/auto_examples/manifold/plot_compare_methods.html#sphx-glr-auto-examples-manifold-plot-compare-methods-py\n# See also:\n# https://scikit-learn.org/stable/auto_examples/manifold/plot_lle_digits.html#\n\n\n\n# To speed-up reduce dimensions by PCA first\n# X_save = X.copy( )\n#r = pca.fit_transform(X)\n#X = r[:1000,:20]\n\n\n\nimport umap \nfrom sklearn import manifold\nfrom sklearn.decomposition import PCA\nfrom sklearn.decomposition import FactorAnalysis\nfrom sklearn.decomposition import NMF\nfrom sklearn.decomposition import FastICA\nfrom sklearn.decomposition import FactorAnalysis\nfrom sklearn.decomposition import LatentDirichletAllocation\nfrom sklearn.ensemble import RandomTreesEmbedding\nfrom sklearn.random_projection import SparseRandomProjection\nfrom sklearn.discriminant_analysis import LinearDiscriminantAnalysis\n\nfrom sklearn.pipeline import make_pipeline\nfrom sklearn.decomposition import TruncatedSVD\n\n\nfrom collections import OrderedDict\nfrom functools import partial\nfrom matplotlib.ticker import NullFormatter\n\n\nn_neighbors = 10\nn_components = 2\n# Set-up manifold methods\nLLE = partial(manifold.LocallyLinearEmbedding,\n              n_neighbors, n_components, eigen_solver='auto')\n\nmethods = OrderedDict()\nmethods['PCA'] = PCA()\nmethods['umap'] = umap.UMAP(n_components = n_components)\nmethods['t-SNE'] = manifold.TSNE(n_components=n_components, init='pca', random_state=0)\nmethods['ICA'] = FastICA(n_components=n_components,         random_state=0)\nmethods['FA'] = FactorAnalysis(n_components=n_components, random_state=0)\n#methods['LLE'] = LLE(method='standard')\n#methods['Modified LLE'] = LLE(method='modified')\n#methods['Isomap'] = manifold.Isomap(n_neighbors, n_components)\nmethods['MDS'] = manifold.MDS(n_components, max_iter=100, n_init=1)\nmethods['SE'] = manifold.SpectralEmbedding(n_components=n_components,\n                                           n_neighbors=n_neighbors)\nmethods['NMF'] = NMF(n_components=n_components,  init='random', random_state=0) \nmethods['RandProj'] = SparseRandomProjection(n_components=n_components, random_state=42)\n\nrand_trees_embed = make_pipeline(RandomTreesEmbedding(n_estimators=200, random_state=0, max_depth=5), TruncatedSVD(n_components=n_components) )\nmethods['RandTrees'] = rand_trees_embed\nmethods['LatDirAll'] = LatentDirichletAllocation(n_components=n_components,  random_state=0)\n#methods['LTSA'] = LLE(method='ltsa') \n#methods['Hessian LLE'] = LLE(method='hessian') \n\nlist_fast_methods = ['PCA','umap','FA', 'NMF','RandProj','RandTrees'] # 'ICA',\nlist_slow_methods = ['t-SNE','LLE','Modified LLE','Isomap','MDS','SE','LatDirAll','LTSA','Hessian LLE']\n\n# transformer = NeighborhoodComponentsAnalysis(init='random',  n_components=2, random_state=0) # Cannot be applied since supervised - requires y \n# methods['LinDisA'] = LinearDiscriminantAnalysis(n_components=n_components)# Cannot be applied since supervised - requires y \n\n\n# Create figure\nfig = plt.figure(figsize=(25, 16))\nplt.suptitle(str_data_inf+ ' n_samples='+str(len(X)),fontsize = 20) \n\n# Plot results\nc = 0\nfor i, (label, method) in enumerate(methods.items()):\n    if label not in list_slow_methods :#  list_fast_methods: #  \n        continue\n        \n    t0 = time.time()\n    try:\n        if isinstance(X, pd.DataFrame):\n            r = method.fit_transform(X.iloc[:N,:])\n        else:\n            r = method.fit_transform(X[:N,:])\n    except:\n        print('Got Exception', label )\n        continue \n    t1 = time.time()\n    print(\"%s: %.2g sec\" % (label, t1 - t0))\n    c+=1\n    fig.add_subplot(2, 3 , c) \n    if isinstance(v4color, pd.Series):\n        sns.scatterplot(x=r[:,0], y=r[:,1] , hue =  v4color.iloc[:N],  palette = palette1, marker =  marker1, alpha = alpha1 )\n    else:\n        sns.scatterplot(x=r[:,0], y=r[:,1] , hue =  v4color[:N],  palette = palette1, marker =  marker1, alpha = alpha1 )\n    plt.title(label,fontsize = 20 )\n    #plt.legend('')\n    if c == 1:\n        plt.setp(ax.get_legend().get_texts(), fontsize=12) # for legend text\n        plt.setp(ax.get_legend().get_title(), fontsize=12) # for legend title    \n        plt.legend( bbox_to_anchor = (1.15, 1.), loc='upper right')\n    else:\n        plt.legend('')\n    \n\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-05-01T15:08:57.615136Z","iopub.execute_input":"2023-05-01T15:08:57.61612Z","iopub.status.idle":"2023-05-01T15:10:49.538143Z","shell.execute_reply.started":"2023-05-01T15:08:57.616079Z","shell.execute_reply":"2023-05-01T15:10:49.53668Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## LocallyLinearEmbedding  (LLE)","metadata":{}},{"cell_type":"code","source":"%%time\n#LLE(n_components=2, method='hessian')\nfrom sklearn.manifold import LocallyLinearEmbedding\n\nfor n_neighbors in [5,10]:\n    print('n_neighbors', n_neighbors,  )\n    # Create figure\n    fig = plt.figure(figsize=(25, 16))\n    plt.suptitle(str_data_inf+ ' n_samples='+str(len(X)),fontsize = 20) \n\n    # Plot results\n    c = 0\n\n\n    for method in ['standard', 'hessian', 'modified', 'ltsa']:\n        if method == 'hessian':\n            # n_components * (n_components + 3) / 2\n            reducer = LocallyLinearEmbedding(n_components=2,  n_neighbors=6,  method = method)\n        else: \n            reducer = LocallyLinearEmbedding(n_components=2, n_neighbors=n_neighbors, method = method)\n\n        t0 = time.time()\n        try:\n            if isinstance(X, pd.DataFrame):\n                r = reducer.fit_transform(X.iloc[:N,:])\n            else:\n                r = reducer.fit_transform(X[:N,:])\n        except Exception as e: \n            print(e)\n            print('Got Exception', method )\n            continue \n        t1 = time.time()\n        print(\"LLE %s: %.2g sec\" % (method, t1 - t0))\n        c+=1\n        fig.add_subplot(2, 2 , c) \n        if isinstance(v4color, pd.Series):\n            sns.scatterplot(x=r[:,0], y=r[:,1] , hue =  v4color.iloc[:N],  palette = palette1, marker =  marker1, alpha = alpha1 )\n        else:\n            sns.scatterplot(x=r[:,0], y=r[:,1] , hue =  v4color[:N],  palette = palette1, marker =  marker1, alpha = alpha1 )\n        plt.title('LLE '+str(method),fontsize = 20 )\n        plt.legend('')\n\n    plt.show() ","metadata":{"execution":{"iopub.status.busy":"2023-05-01T15:10:49.539876Z","iopub.execute_input":"2023-05-01T15:10:49.540615Z","iopub.status.idle":"2023-05-01T15:11:10.360218Z","shell.execute_reply.started":"2023-05-01T15:10:49.540572Z","shell.execute_reply":"2023-05-01T15:11:10.359176Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Isomap","metadata":{}},{"cell_type":"code","source":"%%time\n\nfrom sklearn.manifold import Isomap\n\nfig = plt.figure(figsize=(25, 8))\nplt.suptitle(str_data_inf+ ' n_samples='+str(len(X)),fontsize = 20) \n\nc = 0\nfor n_neighbors in [5,10]:\n    reducer = Isomap(n_components=2,  n_neighbors= n_neighbors)\n\n    t0 = time.time()\n    try:\n        if isinstance(X, pd.DataFrame):\n            r = reducer.fit_transform(X.iloc[:N,:])\n        else:\n            r = reducer.fit_transform(X[:N,:])\n    except Exception as e: \n        print(e)\n        print('Got Exception', n_neighbors )\n        continue \n    t1 = time.time()\n    print(\"Isomap n_neighbors %s: %.2g sec\" % (n_neighbors, t1 - t0))\n    c+=1\n    fig.add_subplot(1, 2 , c) \n    if isinstance(v4color, pd.Series):\n        sns.scatterplot(x=r[:,0], y=r[:,1] , hue =  v4color.iloc[:N],  palette = palette1, marker =  marker1, alpha = alpha1 )\n    else:\n        sns.scatterplot(x=r[:,0], y=r[:,1] , hue =  v4color[:N],  palette = palette1, marker =  marker1, alpha = alpha1 )\n    plt.title('Isomap n_neighbors'+str(n_neighbors),fontsize = 20 )\n    if c == 1:\n        plt.setp(ax.get_legend().get_texts(), fontsize=12) # for legend text\n        plt.setp(ax.get_legend().get_title(), fontsize=12) # for legend title    \n        plt.legend( bbox_to_anchor = (1.15, 1.), loc='upper right')\n    else:\n        plt.legend('')\n\n    \nplt.show() ","metadata":{"execution":{"iopub.status.busy":"2023-05-01T15:11:10.363835Z","iopub.execute_input":"2023-05-01T15:11:10.364617Z","iopub.status.idle":"2023-05-01T15:11:16.037377Z","shell.execute_reply.started":"2023-05-01T15:11:10.36456Z","shell.execute_reply":"2023-05-01T15:11:16.035992Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print('%.1f seconds passed total '%(time.time()-t0start) )","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}