{"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":"# Rudy's TPS July 2022\n\nHere I'm going to attempt a Gaussian Mixture model to cluster the points in our data set.\n\nI haven't tried something like this before, but it's a good time to start!","metadata":{}},{"cell_type":"markdown","source":"**Load the Data**","metadata":{}},{"cell_type":"code","source":"import numpy as np\nimport pandas as pd\nfrom sklearn.preprocessing import PowerTransformer\nfrom sklearn.mixture import BayesianGaussianMixture\nimport seaborn as sns\nimport matplotlib.pyplot as plt\n\nsns.set(rc={'figure.figsize':(10,7)})\n\nsample_submission = pd.read_csv(\"/kaggle/input/tabular-playground-series-jul-2022/sample_submission.csv\")\ndata = pd.read_csv(\"/kaggle/input/tabular-playground-series-jul-2022/data.csv\")\ndata = data.drop(\"id\",axis=1)","metadata":{"execution":{"iopub.status.busy":"2022-07-14T08:40:10.753983Z","iopub.execute_input":"2022-07-14T08:40:10.754692Z","iopub.status.idle":"2022-07-14T08:40:12.411112Z","shell.execute_reply.started":"2022-07-14T08:40:10.754654Z","shell.execute_reply":"2022-07-14T08:40:12.409203Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**Begin with checking out some characteristics of the data set**","metadata":{}},{"cell_type":"code","source":"data","metadata":{"execution":{"iopub.status.busy":"2022-07-14T08:40:12.413612Z","iopub.execute_input":"2022-07-14T08:40:12.414492Z","iopub.status.idle":"2022-07-14T08:40:12.480407Z","shell.execute_reply.started":"2022-07-14T08:40:12.414424Z","shell.execute_reply":"2022-07-14T08:40:12.478941Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"data.info()","metadata":{"execution":{"iopub.status.busy":"2022-07-14T08:40:12.482819Z","iopub.execute_input":"2022-07-14T08:40:12.483316Z","iopub.status.idle":"2022-07-14T08:40:12.512281Z","shell.execute_reply.started":"2022-07-14T08:40:12.483276Z","shell.execute_reply":"2022-07-14T08:40:12.510759Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"*No null values. 28 dimensional data set with 98000 rows.*\n\n*Let's view a correlation map and some distributions for our columns.*","metadata":{}},{"cell_type":"code","source":"sns.heatmap(data.corr(), cmap='Reds')","metadata":{"execution":{"iopub.status.busy":"2022-07-14T08:40:12.513893Z","iopub.execute_input":"2022-07-14T08:40:12.514692Z","iopub.status.idle":"2022-07-14T08:40:13.567469Z","shell.execute_reply.started":"2022-07-14T08:40:12.514634Z","shell.execute_reply":"2022-07-14T08:40:13.566275Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"figure = plt.figure(figsize=(16, 18))\nfeatCount = 29\nfor i in range(featCount):\n    if i < 10:\n        feat_name = 'f_0' + str(i)\n    else:\n        feat_name = 'f_' + str(i)\n    plt.subplot(10, 3, i+1)\n    plt.hist(data[feat_name], bins=100)\n    plt.title(f'{feat_name}')\nfigure.tight_layout(h_pad=1.0, w_pad=1.0)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-07-14T08:40:13.568823Z","iopub.execute_input":"2022-07-14T08:40:13.569196Z","iopub.status.idle":"2022-07-14T08:40:23.412854Z","shell.execute_reply.started":"2022-07-14T08:40:13.569162Z","shell.execute_reply":"2022-07-14T08:40:23.411336Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"*Features 7-13 (our int value columns) have a positive skew distribution. the rest appear to follow a bell curve with mean ~0*\n\nMost of the features don't appear to be very useful for the purposes of clustering. Let's extract only a subset of the features.","metadata":{}},{"cell_type":"markdown","source":"**Scaling our Data**","metadata":{}},{"cell_type":"markdown","source":"Before attempting to fit our model, we need to scale the data. We will use PowerTransform here","metadata":{}},{"cell_type":"markdown","source":"Take only the subset we're interested in","metadata":{}},{"cell_type":"code","source":"cols = ['f_07', 'f_08', 'f_09', 'f_10', 'f_11', 'f_12', 'f_13', 'f_22', 'f_23', 'f_24', 'f_25', 'f_26', 'f_27', 'f_28']","metadata":{"execution":{"iopub.status.busy":"2022-07-14T08:40:23.417351Z","iopub.execute_input":"2022-07-14T08:40:23.417747Z","iopub.status.idle":"2022-07-14T08:40:23.425704Z","shell.execute_reply.started":"2022-07-14T08:40:23.417712Z","shell.execute_reply":"2022-07-14T08:40:23.423522Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"apply the scaler:","metadata":{}},{"cell_type":"code","source":"from sklearn.preprocessing import PowerTransformer\n\ntransformer = PowerTransformer()\nX_scaled = transformer.fit_transform(data[cols])\n\nScaled_Dataframe = pd.DataFrame(X_scaled)\nScaled_Dataframe.describe().T","metadata":{"execution":{"iopub.status.busy":"2022-07-14T08:40:23.427392Z","iopub.execute_input":"2022-07-14T08:40:23.428369Z","iopub.status.idle":"2022-07-14T08:40:25.341850Z","shell.execute_reply.started":"2022-07-14T08:40:23.428324Z","shell.execute_reply":"2022-07-14T08:40:25.340646Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**Principal Component Analysis**","metadata":{}},{"cell_type":"markdown","source":"Conducting PCA on the data will cut down on computational expense by eliminating some unnecessary features from training.\n\nWe're going to reduce our data down to a 3 dimensional data set, while trying to keep as much useful information as possible.\n\nNOTE: I have not used the PCA in the final prediction. It compresses the data too much. I have used just the scaled data of the selected features from above. I have kept the PCA in the notebook as a visualization technique.","metadata":{}},{"cell_type":"code","source":"from sklearn.decomposition import PCA\n\npca = PCA(n_components=3,random_state=1)\n\npca.fit(X_scaled)\n\nX_scaled_pca = pd.DataFrame(pca.transform(X_scaled), columns=([\"col1\",\"col2\",\"col3\"]))\n\nX_scaled_pca.describe().T","metadata":{"_kg_hide-output":true,"execution":{"iopub.status.busy":"2022-07-14T08:40:25.343571Z","iopub.execute_input":"2022-07-14T08:40:25.344297Z","iopub.status.idle":"2022-07-14T08:40:25.981342Z","shell.execute_reply.started":"2022-07-14T08:40:25.344250Z","shell.execute_reply":"2022-07-14T08:40:25.979780Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**Gaussian Mixture Model**","metadata":{}},{"cell_type":"markdown","source":"*Try a gaussian mixture model on our data set. First we will need to find the ideal amount of clusters.*\n\n*Try to find the model that minimizes a theoretical information criterion such as the Bayesian Information Criterion (BIC) or the Akaike Information Criterion (AIC)*","metadata":{}},{"cell_type":"code","source":"from sklearn.mixture import GaussianMixture","metadata":{"execution":{"iopub.status.busy":"2022-07-14T08:40:25.988688Z","iopub.execute_input":"2022-07-14T08:40:25.989795Z","iopub.status.idle":"2022-07-14T08:40:25.997925Z","shell.execute_reply.started":"2022-07-14T08:40:25.989717Z","shell.execute_reply":"2022-07-14T08:40:25.995956Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Try fitting once with 10 components and check the AIC/BIC.","metadata":{}},{"cell_type":"code","source":"gm = GaussianMixture(n_components=10, n_init=10)\ngm.fit(X_scaled_pca)","metadata":{"execution":{"iopub.status.busy":"2022-07-14T08:40:26.000325Z","iopub.execute_input":"2022-07-14T08:40:26.001391Z","iopub.status.idle":"2022-07-14T08:40:54.302292Z","shell.execute_reply.started":"2022-07-14T08:40:26.001322Z","shell.execute_reply":"2022-07-14T08:40:54.301055Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print('Iterations needed: ', gm.n_iter_)","metadata":{"execution":{"iopub.status.busy":"2022-07-14T08:40:54.304307Z","iopub.execute_input":"2022-07-14T08:40:54.305122Z","iopub.status.idle":"2022-07-14T08:40:54.311699Z","shell.execute_reply.started":"2022-07-14T08:40:54.305073Z","shell.execute_reply":"2022-07-14T08:40:54.310517Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print('The AIC score for 10 clusters on this data set is: ', gm.aic(X_scaled_pca))\nprint('The BIC score for 10 clusters on this data set is: ', gm.bic(X_scaled_pca))","metadata":{"execution":{"iopub.status.busy":"2022-07-14T08:40:54.313707Z","iopub.execute_input":"2022-07-14T08:40:54.314540Z","iopub.status.idle":"2022-07-14T08:40:54.586605Z","shell.execute_reply.started":"2022-07-14T08:40:54.314492Z","shell.execute_reply":"2022-07-14T08:40:54.585364Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"We could use this model with 10 clusters to predict the clusters for each of our data points. However it would be preferable to find the optimal number of clusters first.\n\nExtend this to a range of k's","metadata":{}},{"cell_type":"code","source":"gms_per_k = [GaussianMixture(n_components=k, n_init=10, random_state=42).fit(X_scaled_pca)\n             for k in range(1,15)]\n\nbics = [model.bic(X_scaled_pca) for model in gms_per_k]\naics = [model.aic(X_scaled_pca) for model in gms_per_k]","metadata":{"execution":{"iopub.status.busy":"2022-07-14T08:40:54.588529Z","iopub.execute_input":"2022-07-14T08:40:54.589430Z","iopub.status.idle":"2022-07-14T08:46:03.211841Z","shell.execute_reply.started":"2022-07-14T08:40:54.589378Z","shell.execute_reply":"2022-07-14T08:46:03.210083Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.figure(figsize=(10, 4))\nplt.plot(range(0,14), bics, \"bo-\", label=\"BIC\")\nplt.plot(range(0,14), aics, \"go--\", label=\"AIC\")\nplt.xlabel(\"$k$\", fontsize=14)\nplt.ylabel(\"Information Criterion\", fontsize=14)\nplt.axis([0, 13, np.min(aics) - 50, np.max(aics) + 50])\nplt.annotate('Minimum',\n             xy=(7, bics[7]),\n             xytext=(0.3, 0.6),\n             textcoords='figure fraction',\n             fontsize=14,\n             arrowprops=dict(facecolor='black', shrink=0.1)\n            )\nplt.legend()\n\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-07-14T09:01:12.922087Z","iopub.execute_input":"2022-07-14T09:01:12.922464Z","iopub.status.idle":"2022-07-14T09:01:13.175580Z","shell.execute_reply.started":"2022-07-14T09:01:12.922433Z","shell.execute_reply":"2022-07-14T09:01:13.174702Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"The minimum seems to be at 7 according to the BIC. So trying a GMM with 7 clusters seems sensible.\n\nThere is one more parameter which we can try to vary to get the optimal solution: covariance type.\n\nThe code below will compute the best covariance type and the best k value","metadata":{}},{"cell_type":"code","source":"min_bic = np.infty\n\nfor k in range(1, 11):\n    for covariance_type in (\"full\", \"tied\", \"spherical\", \"diag\"):\n        bic = GaussianMixture(n_components=k, n_init=10,\n                              covariance_type=covariance_type,\n                              random_state=42).fit(X_scaled_pca).bic(X_scaled_pca)\n        if bic < min_bic:\n            min_bic = bic\n            best_k = k\n            best_covariance_type = covariance_type","metadata":{"execution":{"iopub.status.busy":"2022-07-14T08:46:03.668181Z","iopub.execute_input":"2022-07-14T08:46:03.668884Z","iopub.status.idle":"2022-07-14T08:53:54.887207Z","shell.execute_reply.started":"2022-07-14T08:46:03.668838Z","shell.execute_reply":"2022-07-14T08:53:54.885718Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print('The optimal number of clusters is: ', best_k)\nprint('The optimal covariance type is: ', best_covariance_type)","metadata":{"execution":{"iopub.status.busy":"2022-07-14T08:53:54.889478Z","iopub.execute_input":"2022-07-14T08:53:54.890321Z","iopub.status.idle":"2022-07-14T08:53:54.903699Z","shell.execute_reply.started":"2022-07-14T08:53:54.890264Z","shell.execute_reply":"2022-07-14T08:53:54.898502Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Let's use these parameters for our GMM","metadata":{}},{"cell_type":"code","source":"gm_pred = GaussianMixture(n_components=7, n_init=10, random_state=42) #'full' is default covariance type\n\n#NOTE: I have reduced the number of components to 7 as after testing, it provides a better fit on the data. \n#Information criterion is a heuristic measure and sometimes doesn't correspond to the best experimental result.\n\npredictions = gm_pred.fit_predict(X_scaled)","metadata":{"execution":{"iopub.status.busy":"2022-07-14T08:53:54.907596Z","iopub.execute_input":"2022-07-14T08:53:54.908486Z","iopub.status.idle":"2022-07-14T08:55:48.728901Z","shell.execute_reply.started":"2022-07-14T08:53:54.908417Z","shell.execute_reply":"2022-07-14T08:55:48.727665Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"We can view a plot of our predicted clusters due to the dimensionality reduction from PCA","metadata":{}},{"cell_type":"code","source":"X_scaled_pca[\"Clusters\"] = predictions\n\nfig = plt.figure(figsize=(10,8))\nax = plt.subplot(111, projection='3d', label=\"bla\")\nax.scatter(X_scaled_pca[\"col1\"], X_scaled_pca[\"col2\"], X_scaled_pca[\"col3\"], s=40, c=X_scaled_pca[\"Clusters\"], marker='o', cmap = 'tab10' )\nax.set_title(\"Cluster plot\")\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-07-14T08:55:48.730520Z","iopub.execute_input":"2022-07-14T08:55:48.730959Z","iopub.status.idle":"2022-07-14T08:55:51.074269Z","shell.execute_reply.started":"2022-07-14T08:55:48.730913Z","shell.execute_reply":"2022-07-14T08:55:51.072215Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"pl = sns.countplot(x=X_scaled_pca[\"Clusters\"])\npl.set_title(\"Distribution Of The Clusters\")\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2022-07-14T08:55:51.075834Z","iopub.execute_input":"2022-07-14T08:55:51.076228Z","iopub.status.idle":"2022-07-14T08:55:51.307978Z","shell.execute_reply.started":"2022-07-14T08:55:51.076192Z","shell.execute_reply":"2022-07-14T08:55:51.307099Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**Bayesian Gaussian Mixture Model**","metadata":{}},{"cell_type":"markdown","source":"Try a new model, this time Bayesian Gaussian\n\nLoad our data as apply transforms as before","metadata":{}},{"cell_type":"code","source":"#sample_submission = pd.read_csv(\"/kaggle/input/tabular-playground-series-jul-2022/sample_submission.csv\")\n#data = pd.read_csv(\"/kaggle/input/tabular-playground-series-jul-2022/data.csv\")\n\n#cols = data.drop(columns=['id']).columns\n#X = data[cols]\n#X_scaled = pd.DataFrame(X, columns=cols)\n\n#scaler = StandardScaler()\n#X_scaled = scaler.fit_transform(X_scaled)\n\n#transformer = PowerTransformer()\n#X_scaled = transformer.fit_transform(X_scaled)\n\n#pca = PCA(n_components=3,random_state=1)\n#pca.fit(X_scaled)\n#X_scaled_pca = pd.DataFrame(pca.transform(X_scaled), columns=([\"col1\",\"col2\",\"col3\"]))\n","metadata":{"execution":{"iopub.status.busy":"2022-07-14T08:55:51.309264Z","iopub.execute_input":"2022-07-14T08:55:51.310155Z","iopub.status.idle":"2022-07-14T08:55:51.315833Z","shell.execute_reply.started":"2022-07-14T08:55:51.310115Z","shell.execute_reply":"2022-07-14T08:55:51.314189Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Now we have reapplied prep steps (with a PowerTransform added in).\n\nTime to try our Bayesian Gaussian Mixture model\n\n","metadata":{}},{"cell_type":"code","source":"#BGM = BayesianGaussianMixture(n_components=7,covariance_type='full',random_state=1,n_init=10)\n\n#predictions = BGM.fit_predict(X_scaled_pca)\n#X_scaled_pca[\"Clusters\"] = predictions\n\n#data = data.drop(\"id\",axis=1)\n#data[\"Clusters\"] = predictions","metadata":{"execution":{"iopub.status.busy":"2022-07-14T08:55:51.317981Z","iopub.execute_input":"2022-07-14T08:55:51.318608Z","iopub.status.idle":"2022-07-14T08:55:51.327990Z","shell.execute_reply.started":"2022-07-14T08:55:51.318556Z","shell.execute_reply":"2022-07-14T08:55:51.327057Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#data","metadata":{"execution":{"iopub.status.busy":"2022-07-14T08:55:51.329502Z","iopub.execute_input":"2022-07-14T08:55:51.330157Z","iopub.status.idle":"2022-07-14T08:55:51.339714Z","shell.execute_reply.started":"2022-07-14T08:55:51.330119Z","shell.execute_reply":"2022-07-14T08:55:51.338202Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#pl = sns.countplot(x=X_scaled_pca[\"Clusters\"])\n#pl.set_title(\"Distribution Of The Clusters\")\n#plt.show()","metadata":{"execution":{"iopub.status.busy":"2022-07-14T08:55:51.341722Z","iopub.execute_input":"2022-07-14T08:55:51.342443Z","iopub.status.idle":"2022-07-14T08:55:51.351575Z","shell.execute_reply.started":"2022-07-14T08:55:51.342381Z","shell.execute_reply":"2022-07-14T08:55:51.350420Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"**Submission**","metadata":{}},{"cell_type":"code","source":"sample_submission.Predicted = pd.DataFrame(predictions)\nsample_submission.to_csv(\"submission.csv\",index=False)","metadata":{"execution":{"iopub.status.busy":"2022-07-14T08:55:51.355399Z","iopub.execute_input":"2022-07-14T08:55:51.356249Z","iopub.status.idle":"2022-07-14T08:55:51.463053Z","shell.execute_reply.started":"2022-07-14T08:55:51.356208Z","shell.execute_reply":"2022-07-14T08:55:51.461994Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}