{"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":"## Clustering of RNA data\n\nPlease visit the article on Medium by Byeungchun Kwon giving all the details: \n\nhttps://medium.com/@byeungchun/machine-learning-for-biological-sequence-data-using-python-573d82f6f17a","metadata":{"id":"hrPaJ3FHSJEB"}},{"cell_type":"markdown","source":" I wondered how much repetition and uneven sampling there is in the training data. Certainly it seems that there is heavy clustering here. ","metadata":{}},{"cell_type":"code","source":"import pandas as pd\nimport os, gc\nimport numpy as np\nfrom sklearn.model_selection import KFold\nimport itertools\n","metadata":{"execution":{"iopub.status.busy":"2023-10-10T23:51:30.042207Z","iopub.execute_input":"2023-10-10T23:51:30.042651Z","iopub.status.idle":"2023-10-10T23:51:32.738244Z","shell.execute_reply.started":"2023-10-10T23:51:30.042607Z","shell.execute_reply":"2023-10-10T23:51:32.737269Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"PATH = '/kaggle/input/stanford-ribonanza-rna-folding-converted/'\nOUT = './'\nbs = 256\nnum_workers = 2\nSEED = 2023\nnfolds = 4\n","metadata":{"execution":{"iopub.status.busy":"2023-10-10T23:51:32.740377Z","iopub.execute_input":"2023-10-10T23:51:32.740845Z","iopub.status.idle":"2023-10-10T23:51:32.746582Z","shell.execute_reply.started":"2023-10-10T23:51:32.740815Z","shell.execute_reply":"2023-10-10T23:51:32.745222Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\nos.makedirs(OUT, exist_ok=True)\ndf = pd.read_parquet(os.path.join(PATH,'train_data.parquet'))","metadata":{"execution":{"iopub.status.busy":"2023-10-10T23:51:32.748457Z","iopub.execute_input":"2023-10-10T23:51:32.748897Z","iopub.status.idle":"2023-10-10T23:51:42.428184Z","shell.execute_reply.started":"2023-10-10T23:51:32.748867Z","shell.execute_reply":"2023-10-10T23:51:42.427308Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"gene_cds = df.head(100000)","metadata":{"execution":{"iopub.status.busy":"2023-10-10T23:51:42.454764Z","iopub.execute_input":"2023-10-10T23:51:42.455544Z","iopub.status.idle":"2023-10-10T23:51:42.460849Z","shell.execute_reply.started":"2023-10-10T23:51:42.455498Z","shell.execute_reply":"2023-10-10T23:51:42.460106Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"##  convert the sequence to k-mer frequency distribution vector","metadata":{"id":"Ik5uOr9UrqSP"}},{"cell_type":"code","source":"from itertools import product","metadata":{"id":"ByiPf9CXrp3P","execution":{"iopub.status.busy":"2023-10-10T23:51:42.461995Z","iopub.execute_input":"2023-10-10T23:51:42.462507Z","iopub.status.idle":"2023-10-10T23:51:42.473099Z","shell.execute_reply.started":"2023-10-10T23:51:42.462478Z","shell.execute_reply":"2023-10-10T23:51:42.472033Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"number_k = 3\nseq_series = list()\nfor id, row in gene_cds.iterrows():\n  seq_id = row['sequence_id']\n  _seq = row['sequence']\n  kmer_tbl = [_seq[i:i+number_k]\n              for i in range(len(_seq)-number_k+1)\n              if _seq[i:i+number_k].find('N') < 0]\n  seq_series.append(\n      pd.Series({x: kmer_tbl.count(x) for x in set(kmer_tbl)}, name=seq_id) / len(kmer_tbl))\n\nseq_series = pd.DataFrame(seq_series).fillna(0.0)","metadata":{"id":"kkGzQ-p-r4qP","execution":{"iopub.status.busy":"2023-10-10T23:51:42.474660Z","iopub.execute_input":"2023-10-10T23:51:42.475005Z","iopub.status.idle":"2023-10-10T23:53:11.412033Z","shell.execute_reply.started":"2023-10-10T23:51:42.474975Z","shell.execute_reply":"2023-10-10T23:53:11.410798Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"seq_series.shape","metadata":{"id":"lLlrZVMED69s","outputId":"ad47e60e-62e6-49ee-97bd-a1f1aec4b028","execution":{"iopub.status.busy":"2023-10-10T23:53:11.414082Z","iopub.execute_input":"2023-10-10T23:53:11.414578Z","iopub.status.idle":"2023-10-10T23:53:11.421883Z","shell.execute_reply.started":"2023-10-10T23:53:11.414523Z","shell.execute_reply":"2023-10-10T23:53:11.421068Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"seq_series.head(3)","metadata":{"id":"j11F4dIfr4hH","outputId":"5ba6cdca-ced7-49b6-f864-2bdd4d17ca59","execution":{"iopub.status.busy":"2023-10-10T23:53:11.423100Z","iopub.execute_input":"2023-10-10T23:53:11.423626Z","iopub.status.idle":"2023-10-10T23:53:11.458579Z","shell.execute_reply.started":"2023-10-10T23:53:11.423596Z","shell.execute_reply":"2023-10-10T23:53:11.457343Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"##  execute k-means clustering","metadata":{"id":"CfusbVJ3Ihwi"}},{"cell_type":"code","source":"import numpy as np\nfrom sklearn.pipeline import Pipeline\nfrom sklearn.preprocessing import LabelEncoder, MinMaxScaler\nfrom sklearn.cluster import KMeans","metadata":{"id":"dYwf9lvOSznf","execution":{"iopub.status.busy":"2023-10-10T23:53:11.461630Z","iopub.execute_input":"2023-10-10T23:53:11.461991Z","iopub.status.idle":"2023-10-10T23:53:11.746426Z","shell.execute_reply.started":"2023-10-10T23:53:11.461955Z","shell.execute_reply":"2023-10-10T23:53:11.745404Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Machine learning pipeline setup\npreprocessor = Pipeline(\n    [\n        (\"scaler\", MinMaxScaler()),\n    ]\n)\nclusterer = Pipeline(\n    [\n        (\n            \"kmeans\",\n            KMeans(\n                n_clusters = 20,\n                init = \"k-means++\",\n                n_init = 50,\n                max_iter = 500,\n                random_state = 42,\n            ),\n        ),\n    ]\n)","metadata":{"id":"Ufh_BBjlx0DC","execution":{"iopub.status.busy":"2023-10-10T23:53:11.748016Z","iopub.execute_input":"2023-10-10T23:53:11.748571Z","iopub.status.idle":"2023-10-10T23:53:11.755245Z","shell.execute_reply.started":"2023-10-10T23:53:11.748525Z","shell.execute_reply":"2023-10-10T23:53:11.754428Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# K-means clustering exercise\nX = seq_series.to_numpy()\npipe = Pipeline(\n    [\n        (\"preprocessor\", preprocessor),\n        (\"clusterer\", clusterer)\n    ]\n)\npipe.fit(X)","metadata":{"id":"OZ8Ntx1mAO6N","outputId":"ce9bfbd9-5ec0-4a47-c784-38a75c883819","execution":{"iopub.status.busy":"2023-10-10T23:53:11.756430Z","iopub.execute_input":"2023-10-10T23:53:11.756728Z","iopub.status.idle":"2023-10-10T23:54:17.602543Z","shell.execute_reply.started":"2023-10-10T23:53:11.756703Z","shell.execute_reply":"2023-10-10T23:54:17.601514Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"_labels = pipe['clusterer']['kmeans'].labels_\nhist = {x: list(_labels).count(x) for x in set(_labels)}\n\nprint(hist)","metadata":{"id":"NG2_Cjg-AOwE","outputId":"83c3dfc0-6da7-4ac6-8d79-26a83aa5a558","execution":{"iopub.status.busy":"2023-10-11T00:05:05.541736Z","iopub.execute_input":"2023-10-11T00:05:05.542139Z","iopub.status.idle":"2023-10-11T00:05:05.707399Z","shell.execute_reply.started":"2023-10-11T00:05:05.542110Z","shell.execute_reply":"2023-10-11T00:05:05.706148Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Visualization","metadata":{"id":"cIol7KjXc6yJ"}},{"cell_type":"code","source":"import matplotlib.pyplot as plt\nimport seaborn as sns\nimport warnings\nwarnings.filterwarnings('ignore')","metadata":{"id":"NVqvjX1mJs5y","execution":{"iopub.status.busy":"2023-10-11T00:05:38.945062Z","iopub.execute_input":"2023-10-11T00:05:38.945506Z","iopub.status.idle":"2023-10-11T00:05:38.952172Z","shell.execute_reply.started":"2023-10-11T00:05:38.945475Z","shell.execute_reply":"2023-10-11T00:05:38.950945Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"fig = plt.figure(figsize = (10, 5))\n\nclusters = list(hist.keys())\ncount = list(hist.values())\n \n# creating the bar plot\nplt.bar(clusters, count, color ='maroon', \n        width = 0.4)\n \nplt.xlabel(\"Cluster\")\nplt.ylabel(\"Count\")\nplt.title(\"Size of clusters\")\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-10-11T00:05:41.865886Z","iopub.execute_input":"2023-10-11T00:05:41.866422Z","iopub.status.idle":"2023-10-11T00:05:42.137170Z","shell.execute_reply.started":"2023-10-11T00:05:41.866387Z","shell.execute_reply":"2023-10-11T00:05:42.135946Z"},"trusted":true},"execution_count":null,"outputs":[]}]}