{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "import seaborn as sns\n",
    "from matplotlib import pyplot as plt\n",
    "\n",
    "from sklearn.preprocessing import StandardScaler\n",
    "from sklearn.model_selection import train_test_split, GridSearchCV, StratifiedKFold, \\\n",
    "learning_curve, validation_curve\n",
    "from sklearn.linear_model import LogisticRegression\n",
    "from sklearn.metrics import accuracy_score, confusion_matrix, precision_score, recall_score\n",
    "from sklearn.decomposition import PCA\n",
    "\n",
    "from scipy.stats import mannwhitneyu\n",
    "\n",
    "sns.set(font_scale=1.5)\n",
    "pd.options.display.max_columns = 50"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Project plan\n",
    "* 1. Feature and data explanation\n",
    "* 2. Primary data analysis\n",
    "* 3. Primary visual data analysis\n",
    "* 4. Insights and found dependencies\n",
    "* 5. Metrics selection\n",
    "* 6. Model selection\n",
    "* 7. Data preprocessing\n",
    "* 8. Cross-validation and adjustment of model hyperparameters\n",
    "* 9. Creation of new features and description of this process\n",
    "* 10. Plotting training and validation curves\n",
    "* 11. Prediction for test or hold-out samples\n",
    "* 12. Conclusions"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. Feature and data explanation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df = pd.read_csv('data.csv')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "df.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1.1 Process of collecting data"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "\n",
    "Features are computed from a digitized image of a fine needle aspirate (FNA) of a breast mass. They describe characteristics of the cell nuclei present in the image. n the 3-dimensional space is that described in: [K. P. Bennett and O. L. Mangasarian: \"Robust Linear Programming Discrimination of Two Linearly Inseparable Sets\", Optimization Methods and Software 1, 1992, 23-34].\n",
    "\n",
    "\n",
    "The data used is available through https://www.kaggle.com/uciml/breast-cancer-wisconsin-data  \n",
    "And can be found on UCI Machine Learning Repository: https://archive.ics.uci.edu/ml/datasets/Breast+Cancer+Wisconsin+%28Diagnostic%29  \n",
    "This database is also available through the UW CS ftp server: ftp ftp.cs.wisc.edu cd math-prog/cpo-dataset/machine-learn/WDBC/"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1.2 Detailed explanation of the task"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The task here is to predict whether the cancer is benign or malignant based in 30 real-valued features."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1.3 Target features"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Attribute Information:\n",
    "\n",
    "1) ID number  \n",
    "2) Diagnosis (M = malignant, B = benign)  \n",
    "3-32) \n",
    "\n",
    "Ten real-valued features are computed for each cell nucleus:  \n",
    "a) radius (mean of distances from center to points on the perimeter)  \n",
    "b) texture (standard deviation of gray-scale values)  \n",
    "c) perimeter  \n",
    "d) area  \n",
    "e) smoothness (local variation in radius lengths)  \n",
    "f) compactness (perimeter^2 / area - 1.0)  \n",
    "g) concavity (severity of concave portions of the contour)  \n",
    "h) concave points (number of concave portions of the contour)  \n",
    "i) symmetry  \n",
    "j) fractal dimension (\"coastline approximation\" - 1)  \n",
    "\n",
    "The mean, standard error and \"worst\" or largest (mean of the three largest values) of these features were computed for each image, resulting in 30 features. For instance, field 3 is Mean Radius, field 13 is Radius SE, field 23 is Worst Radius.\n",
    "\n",
    "\n",
    "All feature values are recoded with four significant digits.  \n",
    "Missing attribute values: none  \n",
    "Class distribution: 357 benign, 212 malignant  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. Primary data analysis"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.0 Data preprocessing"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "target = pd.DataFrame(df['diagnosis'])\n",
    "data = df.drop(['diagnosis'], axis=1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.1 Constant columns"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "General data overview:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "data.info()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Drop constant column ** Unnamed: 32 ** and **id** column which is useless for analize:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "data.drop(['Unnamed: 32', 'id'], axis=1, inplace=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.2 Missing values"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Check data for missing values:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(\"Are there missing values:\", data.isnull().values.any())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.3 Summary statistics"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "General data statistics overview:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "data.describe()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Conclusion:** here we can see vary different min/max values for features, for example *area_mean* and *smoothness_mean*. Thus we should check for outliers (box plot is good option for that)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.4 Statistics for different classes"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Check if difference of features mean values is statistically important. We will use Mann Whitney criteria, because of it unsesetive to outliers and different samples distribution."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for column in data.columns:\n",
    "    m = data[column][target['diagnosis']=='M']\n",
    "    b = data[column][target['diagnosis']=='B']\n",
    "    statistic, pvalue = mannwhitneyu(m, b)\n",
    "    \n",
    "    print('Column:', column, 'Important:', pvalue < 0.05 )"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Conclusion:** differences in almost all features are statistically important. So they will contribute more enough information to classification."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.5 Target feature"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Number of eamples for each class:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "target['diagnosis'].value_counts()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's check the ratio of examples belong to each class:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "target['diagnosis'].value_counts() / target['diagnosis'].size"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Conclusion:** there are a lot more examples for benign class, but not enough for skewed classes problem."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3. Primary visual data analysis "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For the sake of informative data visualization we need to standardize and scale features, because of some features have very different max/min values."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "scaler = StandardScaler()\n",
    "scaled_data = scaler.fit_transform(data)\n",
    "data_scaled = pd.DataFrame(scaled_data, columns=data.columns)\n",
    "data_scaled['diagnosis'] = target['diagnosis']"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 3.1 Linear dependecies of the features (correlation matrix):"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Helper function for plotting feature correlations:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_corr(data):\n",
    "    plt.figure(figsize=[40, 40])\n",
    "    ax = sns.heatmap(data.corr(), annot=True, fmt= '.1f', linewidths=.5)\n",
    "    ax.set_xticklabels(ax.get_xticklabels(), size='xx-large')\n",
    "    ax.set_yticklabels(ax.get_yticklabels(), size='xx-large')\n",
    "    plt.show();"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Data correlations:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_corr(data)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Conclusion:** there are several groups of correlated features:\n",
    "- radius_mean, perimeter_mean, area_mean \n",
    "- compactness_mean, concavity_mean, concave points_mean\n",
    "- radius_se, perimeter_se, area_se\n",
    "- radius_worst, perimeter_worst and area_worst \n",
    "- compactness_worst, concavity_worst, concave points_worst\n",
    "- compactness_se, concavity_se, concave points_se\n",
    "- texture_mean, texture_worst\n",
    "- area_worst, area_mean"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 3.2 Outliers"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "data_z = pd.melt(data_scaled, id_vars=\"diagnosis\", var_name=\"features\", value_name='value')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.figure(figsize=(20, 10));\n",
    "ax = sns.boxplot(x='features', y='value', hue='diagnosis', data=data_z);\n",
    "ax.set_xticklabels(ax.get_xticklabels());\n",
    "plt.xticks(rotation=90);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Conclusion:** there are a lot of variable with outliers. So before training we have to handle it. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 3.3 Distribution of classes"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.figure(figsize=(30, 20));\n",
    "ax = sns.violinplot(x=\"features\", y=\"value\", hue=\"diagnosis\", data=data_z, split=True, inner=\"quartile\");\n",
    "ax.set_xticklabels(ax.get_xticklabels(), size='large');\n",
    "plt.xticks(rotation=90);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Conclusion:** in some features, like *radius_mean*, *texture_mean*, median of each class separated, so they can be useful for classification. Other features, like *smoothness_se*, are not so separated and my be less useful for classification. Most all the features have normal-like distribution with long tail."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 3.4 Dimensionality reduction"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Apply pca for dimensionality reduction:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "pca = PCA(random_state=24)\n",
    "pca.fit(scaled_data)\n",
    "\n",
    "plt.figure(figsize=(10, 10))\n",
    "plt.plot(pca.explained_variance_ratio_, linewidth=2)\n",
    "plt.xlabel('Number of components');\n",
    "plt.ylabel('Explained variance ratio');"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Conclusion:** according to elbow method 3 components may be choosen."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Check the number of components for explaining data variance:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "components = range(1, pca.n_components_ + 1)\n",
    "plt.figure(figsize=(15, 5));\n",
    "plt.bar(components, np.cumsum(pca.explained_variance_ratio_));\n",
    "plt.hlines(y = .95, xmin=0, xmax=len(components), colors='green');"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Conclusion:** The two first components explains the 0.6324 of the variance. We need 10 principal components to explain more than 0.95 of the variance and 17 to explain more than 0.99. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Reduce dimensions of data and plot it:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "pca_two_comp = PCA(n_components=2, random_state=24)\n",
    "two_comp_data = pca_two_comp.fit_transform(scaled_data)\n",
    "plt.scatter(x=two_comp_data[:, 0], y=two_comp_data[:, 1], \n",
    "            c=target['diagnosis'].map({'M': 'red', 'B': 'green'}))\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Conclusion:** data is good enough separable using only two components."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4. Insights and found dependencies"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Data summary:\n",
    "- there are a lot of groups with correlated features. Next we have to get rid from multi-collinearity by selectig one feature for each group.\n",
    "- ration of examples in each class 0.67/0.27. No skewed classes here, which is important for metric selection;\n",
    "- differences in features stitistics (mean) for each class are statisticalli important. So this features will be important for classification.\n",
    "- there are outliers in data. It's important to get rid of them for outliers-sensetive models (logistic regression  for example) before training;\n",
    "- PCA shows thad data is good enough separable using only 3-5 features."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 5. Metrics selection"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Predict whether the cancer is benign or malignant is a **binary classification** task. Here we don't face the probem of skewed classes. So **accuracy** metric will be a good choice for model evaluation. Also this metric is simple enough, thus highly interpretable."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$Accuracy=\\frac{Number~of~corrected~predictions}{Total~number~of~predictions}$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Aslo for the test set we will calculate **precision** and **recall**."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 6. Model selection"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As model was selected **Logistic regression** because:\n",
    "- works well with non categorical features (in our data all features are continious);\n",
    "- robust to small noise in the data;\n",
    "- cases of multi-collinearity can be handled by implementing regularization;\n",
    "- works well if there are no missing data;\n",
    "- efficient implementation availavle;\n",
    "- feature space of current task is not large."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 7. Data preprocessing"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 7.1 Drop useless columns"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Drop constant column **Unnamed: 32** and useless folumn **id** for classification."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X = df.drop(['id', 'Unnamed: 32', 'diagnosis'], axis=1)\n",
    "y = df['diagnosis'].map(lambda x: 1 if x=='M' else 0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 7.3 Split data into train/test\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Split data into train/test with proportional 0.7/0.3 which is common split for such amount of data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.3, random_state=24)\n",
    "print('Train size:', X_train.size)\n",
    "print('Test size:', X_test.size)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 7.2 Feature selection"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "First of all we should handle multi-collinearity. From each group of correleted features we will select only by one feature. So here columns to drop:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "corr_columns = ['perimeter_mean','radius_mean','compactness_mean',\n",
    "                'concave points_mean','radius_se','perimeter_se',\n",
    "                'radius_worst','perimeter_worst','compactness_worst',\n",
    "                'concave points_worst','compactness_se','concave points_se',\n",
    "                'texture_worst','area_worst',\n",
    "                'concavity_mean']"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Drop correlated columns from train data:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X_train = X_train.drop(corr_columns, axis=1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Drop correlated columns from test data:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X_test = X_test.drop(corr_columns, axis=1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Check number of features left:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('Current number of features:', X_train.shape[1])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 7.3 Feature scaling"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "scaler = StandardScaler()\n",
    "X_train_scaled = scaler.fit_transform(X_train)\n",
    "X_test_scaled = scaler.transform(X_test)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 8. Cross-validation and adjustment of model hyperparameters"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Use 3 splits because of we don't have large amount of training data and shuffle samples in random order."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cv = StratifiedKFold(n_splits=3, random_state=24)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Model:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "model = LogisticRegression(random_state=24)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Model parameters:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "model_parameters = {'penalty': ['l1', 'l2'],\n",
    "                    'C': np.linspace(.1, 1, 10)}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To find best hyperparameters we will use grid search as in is quite simple and efficient enough."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "grig_search = GridSearchCV(model, model_parameters, n_jobs=-1, cv=cv, scoring='accuracy')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%time\n",
    "grig_search.fit(X_train_scaled, y_train);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Best model parameters:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "grig_search.best_params_"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Best cv score:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('Accuracy:', grig_search.best_score_)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 9. Creation of new features"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Helper function for applying map operation to data frame attributes:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def apply_cat_op(data, attrs, operation, prefix):\n",
    "    \"\"\"\n",
    "    Apply one operation to data attributes.\n",
    "    \"\"\"\n",
    "    series = [data[attr].map(operation) for attr in attrs]\n",
    "    \n",
    "    _data = pd.concat(series, axis=1).add_prefix(prefix)\n",
    "    new_attrs = _data.columns.values\n",
    "    \n",
    "    return _data, new_attrs"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Creating new features based on medicine requires strong domain knowledge. So we will create them based on mathematics nature of current features. Basic approach for numerical features for regression model is to calculate squares of features in order to capture non-linear dependencies."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Square function:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sq_operation = lambda x: x**2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Create squared feature for each columns and test in with model:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for column in X_train.columns:\n",
    "    X_train_sq, sq_attr = apply_cat_op(X_train, [column], sq_operation, 'sq_')\n",
    "    data = pd.concat([X_train, X_train_sq], axis=1)\n",
    "    \n",
    "    scaler = StandardScaler()\n",
    "    data_scaled = scaler.fit_transform(data)\n",
    "    \n",
    "    grig_search = GridSearchCV(model, model_parameters, n_jobs=-1, cv=cv, scoring='accuracy')\n",
    "    \n",
    "    grig_search.fit(data_scaled, y_train);\n",
    "    \n",
    "    print('Column:', column, ' ', \n",
    "          'Accuracy:', grig_search.best_score_, ' ',\n",
    "          'Best params:', grig_search.best_params_)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As we ca see squaring feature *fractal_dimension_mean*, gives score improving with params {'C': 0.2, 'penalty': 'l2'}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Add new feature to train data:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X_train_sq, atr = apply_cat_op(X_train, ['fractal_dimension_mean'], sq_operation, 'sq_')\n",
    "X_train = pd.concat([X_train, X_train_sq], axis=1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Add new feature to test data:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X_test_sq, atr = apply_cat_op(X_test, ['fractal_dimension_mean'], sq_operation, 'sq_')\n",
    "X_test = pd.concat([X_test, X_test_sq], axis=1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Scale the final data:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "scaler = StandardScaler()\n",
    "X_train_scaled = scaler.fit_transform(X_train)\n",
    "X_test_scaled = scaler.transform(X_test)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Train model  with best parameters on all train data:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "final_model = LogisticRegression(penalty='l2', C=0.2)\n",
    "final_model.fit(X_train_scaled, y_train)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 10. Plotting training and validation curves"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 10.1 Training curve"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Plotting [learning curve fuction](https://scikit-learn.org/stable/auto_examples/model_selection/plot_learning_curve.html#sphx-glr-auto-examples-model-selection-plot-learning-curve-py):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_learning_curve(estimator, title, X, y, ylim=None, cv=None,\n",
    "                        n_jobs=None, train_sizes=np.linspace(.1, 1.0, 5)):\n",
    "    \"\"\"\n",
    "    Generate a simple plot of the test and training learning curve.\n",
    "\n",
    "    Parameters\n",
    "    ----------\n",
    "    estimator : object type that implements the \"fit\" and \"predict\" methods\n",
    "        An object of that type which is cloned for each validation.\n",
    "\n",
    "    title : string\n",
    "        Title for the chart.\n",
    "\n",
    "    X : array-like, shape (n_samples, n_features)\n",
    "        Training vector, where n_samples is the number of samples and\n",
    "        n_features is the number of features.\n",
    "\n",
    "    y : array-like, shape (n_samples) or (n_samples, n_features), optional\n",
    "        Target relative to X for classification or regression;\n",
    "        None for unsupervised learning.\n",
    "\n",
    "    ylim : tuple, shape (ymin, ymax), optional\n",
    "        Defines minimum and maximum yvalues plotted.\n",
    "\n",
    "    cv : int, cross-validation generator or an iterable, optional\n",
    "        Determines the cross-validation splitting strategy.\n",
    "        Possible inputs for cv are:\n",
    "          - None, to use the default 3-fold cross-validation,\n",
    "          - integer, to specify the number of folds.\n",
    "          - :term:`CV splitter`,\n",
    "          - An iterable yielding (train, test) splits as arrays of indices.\n",
    "\n",
    "        For integer/None inputs, if ``y`` is binary or multiclass,\n",
    "        :class:`StratifiedKFold` used. If the estimator is not a classifier\n",
    "        or if ``y`` is neither binary nor multiclass, :class:`KFold` is used.\n",
    "\n",
    "        Refer :ref:`User Guide <cross_validation>` for the various\n",
    "        cross-validators that can be used here.\n",
    "\n",
    "    n_jobs : int or None, optional (default=None)\n",
    "        Number of jobs to run in parallel.\n",
    "        ``None`` means 1 unless in a :obj:`joblib.parallel_backend` context.\n",
    "        ``-1`` means using all processors. See :term:`Glossary <n_jobs>`\n",
    "        for more details.\n",
    "\n",
    "    train_sizes : array-like, shape (n_ticks,), dtype float or int\n",
    "        Relative or absolute numbers of training examples that will be used to\n",
    "        generate the learning curve. If the dtype is float, it is regarded as a\n",
    "        fraction of the maximum size of the training set (that is determined\n",
    "        by the selected validation method), i.e. it has to be within (0, 1].\n",
    "        Otherwise it is interpreted as absolute sizes of the training sets.\n",
    "        Note that for classification the number of samples usually have to\n",
    "        be big enough to contain at least one sample from each class.\n",
    "        (default: np.linspace(0.1, 1.0, 5))\n",
    "    \"\"\"\n",
    "    plt.figure()\n",
    "    plt.title(title)\n",
    "    if ylim is not None:\n",
    "        plt.ylim(*ylim)\n",
    "    plt.xlabel(\"Training examples\")\n",
    "    plt.ylabel(\"Score\")\n",
    "    train_sizes, train_scores, test_scores = learning_curve(\n",
    "        estimator, X, y, cv=cv, n_jobs=n_jobs, train_sizes=train_sizes)\n",
    "    train_scores_mean = np.mean(train_scores, axis=1)\n",
    "    train_scores_std = np.std(train_scores, axis=1)\n",
    "    test_scores_mean = np.mean(test_scores, axis=1)\n",
    "    test_scores_std = np.std(test_scores, axis=1)\n",
    "    plt.grid()\n",
    "\n",
    "    plt.fill_between(train_sizes, train_scores_mean - train_scores_std,\n",
    "                     train_scores_mean + train_scores_std, alpha=0.1,\n",
    "                     color=\"r\")\n",
    "    plt.fill_between(train_sizes, test_scores_mean - test_scores_std,\n",
    "                     test_scores_mean + test_scores_std, alpha=0.1, color=\"g\")\n",
    "    plt.plot(train_sizes, train_scores_mean, 'o-', color=\"r\",\n",
    "             label=\"Training score\")\n",
    "    plt.plot(train_sizes, test_scores_mean, 'o-', color=\"g\",\n",
    "             label=\"Cross-validation score\")\n",
    "\n",
    "    plt.legend(loc=\"best\")\n",
    "    return plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_learning_curve(final_model, 'Logistic regression', \n",
    "                    X_train_scaled, y_train, cv=cv);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Conclusion:** such gap between training and validating curve indicates overfitting. But we can see that validation curve increasing with increasing amount of training examples, so more data is likely to help beat overfitting."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 10.2 Validation curve"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Plotting validation curve function:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_validation_curve(estimator, title, X, y, param_name, param_range, \n",
    "                          cv=None, scoring=None, ylim=None, n_jobs=None):\n",
    "    \"\"\"\n",
    "    Generates a simple plot of training and validation scores for different parameter values.\n",
    "    \n",
    "    Parameters\n",
    "    ----------\n",
    "    estimator : object type that implements the \"fit\" and \"predict\" methods\n",
    "        An object of that type which is cloned for each validation.\n",
    "\n",
    "    title : string\n",
    "        Title for the chart.\n",
    "\n",
    "    X : array-like, shape (n_samples, n_features)\n",
    "        Training vector, where n_samples is the number of samples and\n",
    "        n_features is the number of features.\n",
    "\n",
    "    y : array-like, shape (n_samples) or (n_samples, n_features), optional\n",
    "        Target relative to X for classification or regression;\n",
    "        None for unsupervised learning.\n",
    "    \n",
    "    param_name : string\n",
    "        Name of the parameter that will be varied.\n",
    "\n",
    "    param_range : array-like, shape (n_values,)\n",
    "        The values of the parameter that will be evaluated.\n",
    "\n",
    "    cv : int, cross-validation generator or an iterable, optional\n",
    "        Determines the cross-validation splitting strategy.\n",
    "        Possible inputs for cv are:\n",
    "          - None, to use the default 3-fold cross-validation,\n",
    "          - integer, to specify the number of folds.\n",
    "          - :term:`CV splitter`,\n",
    "          - An iterable yielding (train, test) splits as arrays of indices.\n",
    "\n",
    "        For integer/None inputs, if ``y`` is binary or multiclass,\n",
    "        :class:`StratifiedKFold` used. If the estimator is not a classifier\n",
    "        or if ``y`` is neither binary nor multiclass, :class:`KFold` is used.\n",
    "\n",
    "        Refer :ref:`User Guide <cross_validation>` for the various\n",
    "        cross-validators that can be used here.\n",
    "    \n",
    "    scoring : string, callable or None, optional, default: None\n",
    "        A string (see model evaluation documentation) or\n",
    "        a scorer callable object / function with signature\n",
    "        ``scorer(estimator, X, y)``.\n",
    "    \n",
    "    ylim : tuple, shape (ymin, ymax), optional\n",
    "        Defines minimum and maximum yvalues plotted.\n",
    "\n",
    "    n_jobs : int or None, optional (default=None)\n",
    "        Number of jobs to run in parallel.\n",
    "        ``None`` means 1 unless in a :obj:`joblib.parallel_backend` context.\n",
    "        ``-1`` means using all processors. See :term:`Glossary <n_jobs>`\n",
    "        for more details.\n",
    "    \n",
    "    \"\"\"\n",
    "    train_scores, test_scores = validation_curve(\n",
    "    estimator, X, y, param_name, param_range,\n",
    "    cv=cv, scoring=scoring, n_jobs=n_jobs)\n",
    "    \n",
    "    train_scores_mean = np.mean(train_scores, axis=1)\n",
    "    train_scores_std = np.std(train_scores, axis=1)\n",
    "    test_scores_mean = np.mean(test_scores, axis=1)\n",
    "    test_scores_std = np.std(test_scores, axis=1)\n",
    "\n",
    "    plt.figure()\n",
    "    plt.grid()\n",
    "    plt.title(title)\n",
    "    plt.xlabel(param_name)\n",
    "    plt.ylabel(\"Score\")\n",
    "    if ylim is not None:\n",
    "        plt.ylim(*ylim)\n",
    "    plt.semilogx(param_range, train_scores_mean, 'o-', label=\"Training score\",\n",
    "                 color=\"darkorange\")\n",
    "    plt.fill_between(param_range, train_scores_mean - train_scores_std,\n",
    "                     train_scores_mean + train_scores_std, alpha=0.2,\n",
    "                     color=\"darkorange\")\n",
    "    plt.semilogx(param_range, test_scores_mean, 'o-', label=\"Cross-validation score\",\n",
    "                 color=\"navy\")\n",
    "    plt.fill_between(param_range, test_scores_mean - test_scores_std,\n",
    "                     test_scores_mean + test_scores_std, alpha=0.2,\n",
    "                     color=\"navy\")\n",
    "    plt.legend(loc=\"best\")\n",
    "    return plt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Plot validation curve for model complexity parameter:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_validation_curve(final_model, 'Logistic regression', X_train_scaled, y_train,\n",
    "                             'C', model_parameters['C'], \n",
    "                              cv=cv, scoring='accuracy');"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Conclusion:** gap between training and validating curve indicates overfitting. The best **C** parameter is 0.2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 11. Prediction for test samples"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Make predictions for test samples:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "test_predictions = final_model.predict(X_test_scaled)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Accuracy score:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('Accuracy test score:', accuracy_score(y_test, test_predictions))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Conclusion:** result on the test samples are comparable to the results on cross-validation, even better. Thus our validation scheme is valid."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Confusion matrix:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "test_confusion_matrix = confusion_matrix(test_predictions, y_test);\n",
    "sns.heatmap(test_confusion_matrix, annot=True, fmt='d');"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "From confusion matrix we can see that we have made a few wrong predicions."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Precision:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('Precision:', precision_score(y_test, test_predictions))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Recall:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('Recall:', recall_score(y_test, test_predictions))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 12. Conclusions"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Although we try simple model, it gives 98% accuracy, 98% precision and 97% recall on the test set. There are several (3-5) most important features for classification, which could indicates that our data is not representable or biased. So, it's a good option to try model on more data. Feature generation based on medicine knowledge for such data is quite challenging, so we build them based on math nature. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Ways of improving:\n",
    "- collect more data and re-train model on it, as we can see validation score improvement with data amount increasing on learning curve;\n",
    "- dig into domain and generate more features based on medicine;\n",
    "- try another models, like neural network (for capturing complex non-linear dependences) or random forest (robust to overfitting);\n",
    "- apply PCA for data dimensionality reduction and train model on reduced data;\n",
    "- try stacking differet models."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}