{"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":"## original notebook link: https://www.kaggle.com/code/yusaku5739/visualization-of-gene-expression-using-umap","metadata":{}},{"cell_type":"markdown","source":"# Gene expression visualization using UMAP\n\n>UMAP (Uniform Manifold Approximation and Projection) is a novel manifold learning technique for dimension reduction. UMAP is constructed from a theoretical framework based in Riemannian geometry and algebraic topology. The result is a practical scalable algorithm that is applicable to real world data. The UMAP algorithm is competitive with t-SNE for visualization quality, and arguably preserves more of the global structure with superior run time performance. Furthermore, UMAP has no computational restrictions on embedding dimension, making it viable as a general purpose dimension reduction technique for machine learning\n(McInnes, Leland; Healy, John; Melville, James (2018-12-07). \"Uniform manifold approximation and projection for dimension reduction\". arXiv:1802.03426.)\n\n### Please Upvote if you Find this Useful :)","metadata":{}},{"cell_type":"code","source":"!pip install umap-learn","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2023-10-02T16:59:50.099295Z","iopub.execute_input":"2023-10-02T16:59:50.099721Z","iopub.status.idle":"2023-10-02T17:00:00.363080Z","shell.execute_reply.started":"2023-10-02T16:59:50.099672Z","shell.execute_reply":"2023-10-02T17:00:00.361782Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\nimport seaborn as sns\nsns.set_theme()\nimport sklearn\nimport umap\nimport random\nfrom sklearn.svm import LinearSVR\nfrom sklearn.linear_model import Ridge\nfrom sklearn.multioutput import MultiOutputRegressor\nfrom sklearn.model_selection import train_test_split","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:00:00.366345Z","iopub.execute_input":"2023-10-02T17:00:00.366944Z","iopub.status.idle":"2023-10-02T17:00:00.376806Z","shell.execute_reply.started":"2023-10-02T17:00:00.366893Z","shell.execute_reply":"2023-10-02T17:00:00.375242Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df = pd.read_parquet(\"/kaggle/input/open-problems-single-cell-perturbations/de_train.parquet\")\ndf_gene_exp = df.drop(columns=[\"cell_type\",\"sm_name\", \"sm_lincs_id\", \"SMILES\", \"control\"])\ndf_info = df[[\"cell_type\",\"sm_name\", \"sm_lincs_id\", \"SMILES\", \"control\"]]","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:00:00.379087Z","iopub.execute_input":"2023-10-02T17:00:00.379730Z","iopub.status.idle":"2023-10-02T17:00:01.797597Z","shell.execute_reply.started":"2023-10-02T17:00:00.379645Z","shell.execute_reply":"2023-10-02T17:00:01.796485Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# do UMAP clustering\nmapper = umap.UMAP(random_state=42,\n                   n_neighbors=5,\n                   min_dist=0.4,\n                   metric=\"correlation\")\nembedding = mapper.fit_transform(df_gene_exp)\n\nembedding_x = embedding[:, 0]\nembedding_y = embedding[:, 1]\n\n\n# plot UMAP\ncolors=[\"r\", \"b\", \"g\", \"y\", \"m\", \"c\", \"k\", \"w\"]\ndic_c = {}\n\nfor (cell,c) in zip(df_info[\"cell_type\"].unique(), colors):\n    dic_c[cell]=c\n    cell_i = df_info[df_info[\"cell_type\"]==cell].index.to_list()\n    plt.scatter(embedding_x[cell_i], embedding_y[cell_i], label=cell, s=10, color=c)\n    \ncont = df_info[df_info[\"control\"]==True]\nfor cell in cont[\"cell_type\"].unique():\n    i = cont.query(\"cell_type==@cell\").index.to_list()\n    plt.scatter(embedding_x[i], embedding_y[i], label=f\"positive_ctrl_{cell}\", marker=\"^\",color=dic_c[cell], s=40, edgecolors=\"black\")\n\nplt.grid()\nplt.legend(loc='upper left', bbox_to_anchor=(1, 1))\n#plt.title()\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:00:01.798948Z","iopub.execute_input":"2023-10-02T17:00:01.799261Z","iopub.status.idle":"2023-10-02T17:00:16.954062Z","shell.execute_reply.started":"2023-10-02T17:00:01.799236Z","shell.execute_reply":"2023-10-02T17:00:16.952803Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_info['cell_type'].unique()","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:00:16.956932Z","iopub.execute_input":"2023-10-02T17:00:16.957841Z","iopub.status.idle":"2023-10-02T17:00:16.965582Z","shell.execute_reply.started":"2023-10-02T17:00:16.957806Z","shell.execute_reply":"2023-10-02T17:00:16.964709Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Explore DE residuals in 'Clotrimazole' group","metadata":{}},{"cell_type":"code","source":"def plot_by_molecule_scatter(molecule, de=df_gene_exp, info=df_info):\n    my_index = info.query(\"sm_name==@molecule\").index\n    subdf_exp = de.loc[de.index.isin(my_index)]\n    subdf_info = info.loc[info.index.isin(my_index)]\n\n    display(subdf_exp.head())\n    display(subdf_info.head())\n    \n    aggr_de = subdf_exp.aggregate([np.mean, np.std])\n    sns.scatterplot(x=aggr_de.loc[\"mean\"], y=aggr_de.loc[\"std\"])\n    plt.show()\n    \n    return aggr_de","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:00:16.966932Z","iopub.execute_input":"2023-10-02T17:00:16.967897Z","iopub.status.idle":"2023-10-02T17:00:16.980559Z","shell.execute_reply.started":"2023-10-02T17:00:16.967841Z","shell.execute_reply":"2023-10-02T17:00:16.979631Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"_ = plot_by_molecule_scatter(\"Clotrimazole\")","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:00:16.982201Z","iopub.execute_input":"2023-10-02T17:00:16.982959Z","iopub.status.idle":"2023-10-02T17:00:29.869574Z","shell.execute_reply.started":"2023-10-02T17:00:16.982905Z","shell.execute_reply":"2023-10-02T17:00:29.868357Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Explore dipersion in cell types","metadata":{}},{"cell_type":"code","source":"def plot_by_cell_type_scatter(cell_type, de=df_gene_exp, info=df_info):\n    print(\"=\"*20)\n    print(cell_type)\n    print(\"=\"*20)\n    \n    my_index = info.query(\"cell_type==@cell_type\").index\n    subdf_exp = de.loc[de.index.isin(my_index)]\n    subdf_info = info.loc[info.index.isin(my_index)]\n\n#     display(subdf_exp.head())\n#     display(subdf_info.head())\n    \n    aggr_de = subdf_exp.aggregate([np.mean, np.std])\n#     sns.scatterplot(x=aggr_de.loc[\"mean\"].abs(), y=aggr_de.loc[\"std\"])\n    sns.regplot(x=aggr_de.loc[\"mean\"].abs(), y=aggr_de.loc[\"std\"],\n                lowess=False,\n                scatter_kws={\"color\": \"blue\"}, line_kws={\"color\": \"red\"})\n    plt.show()\n    \n    return aggr_de","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:00:29.871295Z","iopub.execute_input":"2023-10-02T17:00:29.873214Z","iopub.status.idle":"2023-10-02T17:00:29.880031Z","shell.execute_reply.started":"2023-10-02T17:00:29.873137Z","shell.execute_reply":"2023-10-02T17:00:29.878724Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for cell_type in df_info['cell_type'].unique():\n    _ = plot_by_cell_type_scatter(cell_type, de=df_gene_exp.iloc[:, :], info=df_info.iloc[:, :])","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:00:29.881476Z","iopub.execute_input":"2023-10-02T17:00:29.881980Z","iopub.status.idle":"2023-10-02T17:01:58.889551Z","shell.execute_reply.started":"2023-10-02T17:00:29.881948Z","shell.execute_reply":"2023-10-02T17:01:58.888645Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## p-value means is correlated with p-value stds across genes in all cell types (but T cells CD8+ scatterplot looks weird). We can calculate overdispersion of p-values to find most overdispersed genes.\n\n* Predict dispersion p-values  for each cell type\n* Find difference between real dispersion and predicted dispersion (overdispersion)\n* Use overdispersion as a feature or filter criterion ","metadata":{}},{"cell_type":"code","source":"from sklearn.linear_model import LinearRegression","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:01:58.890968Z","iopub.execute_input":"2023-10-02T17:01:58.891558Z","iopub.status.idle":"2023-10-02T17:01:58.895916Z","shell.execute_reply.started":"2023-10-02T17:01:58.891525Z","shell.execute_reply":"2023-10-02T17:01:58.894971Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def find_overdispersion(cell_type, de=df_gene_exp, info=df_info):\n    \n    # subset data\n    my_index = info.query(\"cell_type==@cell_type\").index\n    subdf_exp = de.loc[de.index.isin(my_index)]\n    subdf_info = info.loc[info.index.isin(my_index)]\n    \n    # find real std\n    std = subdf_exp.aggregate(np.std)\n    \n    # prepare feature df (transpose to sampleXfeature matrix and abs() pvals)\n    features = subdf_exp.abs().T\n    \n    # fit linear regression\n    reg = LinearRegression().fit(features, std)\n    print(\"=\"*30)\n    print(f\"Linear regression score for {cell_type}: {reg.score(features, std)}\")\n    print(\"=\"*30)\n    # predict std\n    predictions = reg.predict(features)\n    \n    # find overdisperssion\n    overdisp = predictions - std\n    \n    \n    \n    fig, axs = plt.subplots(ncols=3)\n    fig.set_size_inches(18.5, 8)\n    sns.regplot(x=subdf_exp.mean(), y=overdisp, ax=axs[0], line_kws={\"color\": \"red\"})\n    sns.histplot(overdisp, ax=axs[1])\n    sns.boxplot(overdisp, ax=axs[2])\n    plt.show()\n    \n    return overdisp","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:01:58.897395Z","iopub.execute_input":"2023-10-02T17:01:58.897839Z","iopub.status.idle":"2023-10-02T17:01:58.910467Z","shell.execute_reply.started":"2023-10-02T17:01:58.897788Z","shell.execute_reply":"2023-10-02T17:01:58.909603Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for cell_type in df_info['cell_type'].unique():\n    locals()[cell_type.replace(\" \", \"_\").replace(\"+\", \"\")] = find_overdispersion(cell_type)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:01:58.911917Z","iopub.execute_input":"2023-10-02T17:01:58.912819Z","iopub.status.idle":"2023-10-02T17:02:16.053391Z","shell.execute_reply.started":"2023-10-02T17:01:58.912786Z","shell.execute_reply":"2023-10-02T17:02:16.051966Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"We got rid of linear dependence of p-values means and pvalue-stds. But there are a lot of overdispersed gene outliers. Lets see if they repeated in each cell type?","metadata":{}},{"cell_type":"code","source":"from matplotlib_venn import venn2","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:02:16.055311Z","iopub.execute_input":"2023-10-02T17:02:16.055762Z","iopub.status.idle":"2023-10-02T17:02:16.076043Z","shell.execute_reply.started":"2023-10-02T17:02:16.055719Z","shell.execute_reply":"2023-10-02T17:02:16.074927Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def get_ranges(gene_list):\n    # finding the 1st quartile\n    q1 = np.quantile(gene_list, 0.25)\n\n    # finding the 3rd quartile\n    q3 = np.quantile(gene_list, 0.75)\n    med = np.median(gene_list)\n\n    # finding the iqr region\n    iqr = q3-q1\n\n    # finding upper and lower whiskers\n    upper_bound = q3+(1.5*iqr)\n    lower_bound = q1-(1.5*iqr)\n    print(f\"iqr:{iqr}, upper_bound:{upper_bound}, lower_bound:{lower_bound}\")\n    \n    outliers = gene_list[(gene_list <= lower_bound) | (gene_list >= upper_bound)]\n    \n    return outliers","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:02:16.080648Z","iopub.execute_input":"2023-10-02T17:02:16.081061Z","iopub.status.idle":"2023-10-02T17:02:16.088744Z","shell.execute_reply.started":"2023-10-02T17:02:16.081031Z","shell.execute_reply":"2023-10-02T17:02:16.086832Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"NK_cells_otliers = get_ranges(NK_cells)\nT_cells_CD4_otliers = get_ranges(T_cells_CD4)\nT_cells_CD8_otliers = get_ranges(T_cells_CD8)\nT_regulatory_cells_otliers = get_ranges(T_regulatory_cells)\nB_cells_otliers = get_ranges(B_cells)\nMyeloid_cells_otliers = get_ranges(Myeloid_cells)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:02:16.090615Z","iopub.execute_input":"2023-10-02T17:02:16.091004Z","iopub.status.idle":"2023-10-02T17:02:16.119429Z","shell.execute_reply.started":"2023-10-02T17:02:16.090975Z","shell.execute_reply":"2023-10-02T17:02:16.118169Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print(\"All intersections\")\nprint(set(NK_cells_otliers.index) & set(T_cells_CD4_otliers.index) & set(T_cells_CD8_otliers.index) & \n      set(T_regulatory_cells_otliers.index)& set(B_cells_otliers.index) & set(Myeloid_cells_otliers.index))","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:02:16.121638Z","iopub.execute_input":"2023-10-02T17:02:16.122164Z","iopub.status.idle":"2023-10-02T17:02:16.130497Z","shell.execute_reply.started":"2023-10-02T17:02:16.122122Z","shell.execute_reply":"2023-10-02T17:02:16.129288Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Intersections in neigbors by UMAP","metadata":{}},{"cell_type":"code","source":"venn2([set(NK_cells_otliers.index), set(T_cells_CD4_otliers.index)], set_labels =[\"NK\", \"CD4\"])","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:02:16.132178Z","iopub.execute_input":"2023-10-02T17:02:16.132620Z","iopub.status.idle":"2023-10-02T17:02:16.272684Z","shell.execute_reply.started":"2023-10-02T17:02:16.132580Z","shell.execute_reply":"2023-10-02T17:02:16.270972Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"venn2([set(T_regulatory_cells_otliers.index), set(T_cells_CD8_otliers.index)], set_labels =[\"T-reg\", \"CD8\"])","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:02:16.274960Z","iopub.execute_input":"2023-10-02T17:02:16.275937Z","iopub.status.idle":"2023-10-02T17:02:16.433083Z","shell.execute_reply.started":"2023-10-02T17:02:16.275852Z","shell.execute_reply":"2023-10-02T17:02:16.431472Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"venn2([set(Myeloid_cells_otliers.index), set(T_cells_CD4_otliers.index)], set_labels =[\"Myeloid\", \"CD4\"])","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:02:16.435222Z","iopub.execute_input":"2023-10-02T17:02:16.436721Z","iopub.status.idle":"2023-10-02T17:02:16.633774Z","shell.execute_reply.started":"2023-10-02T17:02:16.436660Z","shell.execute_reply":"2023-10-02T17:02:16.632088Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Looks like there are independent genes with high dispersion in each cell type","metadata":{}},{"cell_type":"code","source":"concat_ods = pd.concat([NK_cells, T_cells_CD4, T_cells_CD8, T_regulatory_cells, B_cells, Myeloid_cells], axis=1).T\nconcat_ods.index = ['NK cells', 'T cells CD4+', 'T cells CD8+', 'T regulatory cells',\n                    'B cells', 'Myeloid cells']\nprint(concat_ods.shape)\nconcat_ods.head(6)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:02:16.636028Z","iopub.execute_input":"2023-10-02T17:02:16.637216Z","iopub.status.idle":"2023-10-02T17:02:16.708065Z","shell.execute_reply.started":"2023-10-02T17:02:16.637147Z","shell.execute_reply":"2023-10-02T17:02:16.707150Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## PCA of overdispersion","metadata":{}},{"cell_type":"code","source":"from sklearn.decomposition import PCA\npca = PCA(n_components=2)\nod_pca = pca.fit_transform(concat_ods)\nod_pca = pd.DataFrame(od_pca, columns = ['PC1', 'PC2'])\n\nplt.figure()\nplt.figure(figsize=(10,10))\nplt.xticks(fontsize=12)\nplt.yticks(fontsize=14)\nplt.xlabel('Principal Component - 1',fontsize=20)\nplt.ylabel('Principal Component - 2',fontsize=20)\ntargets = [0,1,2,3,4,5]\ncolors = [\"r\", \"b\", \"g\", \"y\", \"m\", \"c\"]\nfor target, color in zip(targets,colors):\n    plt.scatter(od_pca.loc[target, 'PC1']\n               , od_pca.loc[target, 'PC2'], c = color, s = 50)\nplt.legend(['NK_cells', 'T_cells_CD4', 'T_cells_CD8', 'T_regulatory_cells', 'B_cells', 'Myeloid_cells'])\nplt.show()\n# NK_cells, T_cells_CD4, T_cells_CD8, T_regulatory_cells, B_cells, Myeloid_cells","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:02:16.709338Z","iopub.execute_input":"2023-10-02T17:02:16.709839Z","iopub.status.idle":"2023-10-02T17:02:17.416075Z","shell.execute_reply.started":"2023-10-02T17:02:16.709809Z","shell.execute_reply":"2023-10-02T17:02:17.414929Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Use overdisperssion to train the model.\n\n__WARNING:__ There is a leak in training data since we count overdispersion on full set of DE data. ","metadata":{}},{"cell_type":"code","source":"id_map = pd.read_csv('../input/open-problems-single-cell-perturbations/id_map.csv')\nid_map.head()","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:02:17.417827Z","iopub.execute_input":"2023-10-02T17:02:17.418649Z","iopub.status.idle":"2023-10-02T17:02:17.442136Z","shell.execute_reply.started":"2023-10-02T17:02:17.418614Z","shell.execute_reply":"2023-10-02T17:02:17.441262Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sample_submission = pd.read_csv('../input/open-problems-single-cell-perturbations/sample_submission.csv')\nprint(sample_submission.shape)\nsample_submission.head()","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:02:17.443216Z","iopub.execute_input":"2023-10-02T17:02:17.443831Z","iopub.status.idle":"2023-10-02T17:02:20.336798Z","shell.execute_reply.started":"2023-10-02T17:02:17.443800Z","shell.execute_reply":"2023-10-02T17:02:20.335797Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Prepare test and train data according to https://www.kaggle.com/code/dangnguyen97/linearsvr","metadata":{}},{"cell_type":"code","source":"xlist  = ['cell_type','sm_name']\n_ylist = ['cell_type','sm_name','sm_lincs_id','SMILES','control']\n\ny = df.drop(columns=_ylist)\ny.shape","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:02:20.338119Z","iopub.execute_input":"2023-10-02T17:02:20.339008Z","iopub.status.idle":"2023-10-02T17:02:20.391327Z","shell.execute_reply.started":"2023-10-02T17:02:20.338975Z","shell.execute_reply":"2023-10-02T17:02:20.390203Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"np.random.seed(42)\ntrain = df[xlist]\ntrain = pd.merge(train, concat_ods, left_on='cell_type', right_index=True)\ntrain = pd.get_dummies(train, columns=xlist, drop_first=True)\nnoise = np.random.normal(0, 0.1, [train.shape[0], concat_ods.shape[1]-1])\ntrain.iloc[:,:concat_ods.shape[1]-1] += noise\ntrain.head()","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:02:20.393812Z","iopub.execute_input":"2023-10-02T17:02:20.394239Z","iopub.status.idle":"2023-10-02T17:02:22.390713Z","shell.execute_reply.started":"2023-10-02T17:02:20.394190Z","shell.execute_reply":"2023-10-02T17:02:22.389110Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"test = pd.merge(id_map[xlist], concat_ods, left_on='cell_type', right_index=True)\ntest = pd.get_dummies(test, columns=xlist, drop_first=True)\nnoise = np.random.normal(0, 0.1, [test.shape[0], concat_ods.shape[1]-1])\ntest.iloc[:,:concat_ods.shape[1]-1] += noise\ntest.head()","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:02:22.392270Z","iopub.execute_input":"2023-10-02T17:02:22.392616Z","iopub.status.idle":"2023-10-02T17:02:24.076187Z","shell.execute_reply.started":"2023-10-02T17:02:22.392589Z","shell.execute_reply":"2023-10-02T17:02:24.074956Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"uncommon = [f for f in train.iloc[:,concat_ods.shape[1]:] if f not in test.iloc[:,concat_ods.shape[1]:]]\n\nX = train.drop(columns=uncommon)\nX.shape","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:02:24.078082Z","iopub.execute_input":"2023-10-02T17:02:24.078832Z","iopub.status.idle":"2023-10-02T17:02:24.154419Z","shell.execute_reply.started":"2023-10-02T17:02:24.078789Z","shell.execute_reply":"2023-10-02T17:02:24.153307Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def mrrmse(y_pred: pd.DataFrame, y_true: pd.DataFrame):\n    \n    return ((y_pred - y_true)**2).mean(axis=1).apply(np.sqrt).mean()","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:02:24.155761Z","iopub.execute_input":"2023-10-02T17:02:24.156126Z","iopub.status.idle":"2023-10-02T17:02:24.162609Z","shell.execute_reply.started":"2023-10-02T17:02:24.156097Z","shell.execute_reply":"2023-10-02T17:02:24.161331Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.20, random_state=421)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:02:24.164712Z","iopub.execute_input":"2023-10-02T17:02:24.165532Z","iopub.status.idle":"2023-10-02T17:02:24.306810Z","shell.execute_reply.started":"2023-10-02T17:02:24.165487Z","shell.execute_reply":"2023-10-02T17:02:24.305657Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\n# model = LinearSVR(max_iter= 1000, epsilon= 0.1)\nmodel = Ridge(random_state=42)\n# wrapper = MultiOutputRegressor(model)\n# wrapper.fit(X_train, y_train)\nmodel.fit(X_train, y_train)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:02:24.308256Z","iopub.execute_input":"2023-10-02T17:02:24.308881Z","iopub.status.idle":"2023-10-02T17:02:32.078363Z","shell.execute_reply.started":"2023-10-02T17:02:24.308827Z","shell.execute_reply":"2023-10-02T17:02:32.076748Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# predict = wrapper.predict(X_test)\npredict = model.predict(X_test)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:02:32.080813Z","iopub.execute_input":"2023-10-02T17:02:32.082950Z","iopub.status.idle":"2023-10-02T17:02:34.517635Z","shell.execute_reply.started":"2023-10-02T17:02:32.082630Z","shell.execute_reply":"2023-10-02T17:02:34.516070Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#mrrmse(pd.DataFrame(predict), pd.DataFrame(y_test.values))\nN = random.randrange(predict.shape[1])\n\nplt.style.use('seaborn-whitegrid') \nplt.figure(figsize=(8, 4), facecolor='lightyellow')\nplt.title(f'Column:  #{N}', fontsize=12)\nplt.gca().set_facecolor('lightgray')\n\nsns.histplot(y_test.values[:,N]-predict[:,N], bins=100, color='red',\n             )\nplt.legend(['y_true','y_pred'], loc=1)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:02:34.524997Z","iopub.execute_input":"2023-10-02T17:02:34.526396Z","iopub.status.idle":"2023-10-02T17:02:35.076499Z","shell.execute_reply.started":"2023-10-02T17:02:34.526314Z","shell.execute_reply":"2023-10-02T17:02:35.074935Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Final model","metadata":{}},{"cell_type":"code","source":"# model = LinearSVR(max_iter= 1500, epsilon= 0.075)\n\n# wrapper = MultiOutputRegressor(model)\nmodel.fit(X, y)","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:02:35.078087Z","iopub.execute_input":"2023-10-02T17:02:35.078592Z","iopub.status.idle":"2023-10-02T17:02:43.807549Z","shell.execute_reply.started":"2023-10-02T17:02:35.078547Z","shell.execute_reply":"2023-10-02T17:02:43.805844Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# submission1 = pd.DataFrame(wrapper.predict(test), columns= df.columns[5:])\nsubmission1 = pd.DataFrame(model.predict(test), columns= df.columns[5:])\nsubmission1.index.name = 'id'\nsubmission1","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:02:43.812793Z","iopub.execute_input":"2023-10-02T17:02:43.817375Z","iopub.status.idle":"2023-10-02T17:02:47.062096Z","shell.execute_reply.started":"2023-10-02T17:02:43.817297Z","shell.execute_reply":"2023-10-02T17:02:47.060262Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"submission1.to_csv('submission1.csv')","metadata":{"execution":{"iopub.status.busy":"2023-10-02T17:02:47.064380Z","iopub.execute_input":"2023-10-02T17:02:47.065503Z","iopub.status.idle":"2023-10-02T17:02:55.902114Z","shell.execute_reply.started":"2023-10-02T17:02:47.065439Z","shell.execute_reply":"2023-10-02T17:02:55.900879Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}