Neural_networks.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "import seaborn as sb\n",
    "import matplotlib.pyplot as plt\n",
    "import pandas as pd\n",
    "import sklearn as sk\n",
    "import numpy as np\n",
    "\n",
    "from sklearn.linear_model import LogisticRegression\n",
    "from sklearn.metrics import mean_squared_error\n",
    "from sklearn import metrics\n",
    "from sklearn.model_selection import cross_validate\n",
    "from sklearn.utils import resample\n",
    "from sklearn.neural_network import MLPClassifier\n",
    "from sklearn.model_selection import train_test_split\n",
    "\n",
    "## read data \n",
    "CIP_data = pd.read_csv(\"CIP_data_encode_prev_not_dropped.csv\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Step 1: create model and calculate apparent performance metric of interest (P)\n",
    "X = CIP_data[['MSMW', 'MSW', 'Oth/Unk/Missing','Northeast', 'Southeast', 'Southwest', 'West', 'PREV_CLINIC', 'PREV_REGION']]\n",
    "y = CIP_data['Susceptible']\n",
    "print(y.sum()/len(y))\n",
    "\n",
    "\n",
    "#X_train, X_test, y_train, y_test = train_test_split(X, y, test_size = 0.30) ## even doing trian/test split doesn't work...\n",
    "\n",
    "neural_network_model = MLPClassifier(solver = 'lbfgs', activation = 'tanh', max_iter = 1000 ,hidden_layer_sizes= 44, random_state=10)\n",
    "\n",
    "neural_network_fit = neural_network_model.fit(X,y)\n",
    "y_predict_nn = neural_network_fit.predict(X)\n",
    "\n",
    "ROC_AUC_neural_network = metrics.roc_auc_score(y, y_predict_nn)\n",
    "\n",
    "print(ROC_AUC_neural_network)\n",
    "\n",
    "print(\"ACCURACY OF THE MODEL: \", metrics.accuracy_score(y, y_predict_nn)) ## but get model accuracy of 0.9470427693866846... so why is the ROC so low?  \n",
    "## any why is this getting the same score as the random forest model?\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "#### now try bootstrapping\n",
    "#bootstrap data\n",
    "n_iterations = 10\n",
    "bootstrapped_stats = []\n",
    "\n",
    "#1. Create model using all data and get ROC_AUC (\"ROC_AUC_neural_network\")\n",
    "X = CIP_data[['MSMW', 'MSW', 'Oth/Unk/Missing','Northeast', 'Southeast', 'Southwest', 'West', 'PREV_REGION', 'PREV_CLINIC', 'DELTA_REGION']]]\n",
    "y = CIP_data['Susceptible']\n",
    "\n",
    "model = MLPClassifier(solver = 'lbfgs', activation = 'tanh', max_iter = 1000 ,hidden_layer_sizes= 4, random_state=10)\n",
    "model_fit = model.fit(X, y)\n",
    "y_predict = model_fit.predict(X)\n",
    "\n",
    "ROC_AUC_neural_network_apparent = metrics.roc_auc_score(y, y_predict)\n",
    "\n",
    "for i in range(n_iterations):\n",
    "       #2. (A) Sample all individuals w/replacement\n",
    "        sample = CIP_data.sample(frac = 1, replace=True) ##(a) sample n individuals with replacement\n",
    "        X_sample = sample[['MSMW', 'MSW', 'Oth/Unk/Missing','Northeast', 'Southeast', 'Southwest', 'West', 'PREV_REGION', 'PREV_CLINIC', 'DELTA_REGION']]]\n",
    "        y_sample = sample['Susceptible']\n",
    "\n",
    "\n",
    "       #  (B) Develop predictive model and find apparent performance\n",
    "        model_fit = model.fit(X_sample, y_sample)\n",
    "        y_sample_predict = model_fit.predict(X_sample)\n",
    "        ROC_AUC_neural_network_bootstrap_sample_performance = metrics.roc_auc_score(y_sample, y_sample_predict) \n",
    "\n",
    "       #  (C) Performance of predictive model on original sample (i.e. original population, X)\n",
    "        y_test_predict = model_fit.predict(X)\n",
    "        ROC_AUC_neural_network_bootstrap_test_performance = metrics.roc_auc_score(y, y_test_predict) ## 0.756384214489288\n",
    "      ### (D) Calculate optimisation by getting (B) - (D) \n",
    "        optimism = ROC_AUC_neural_network_bootstrap_sample_performance - ROC_AUC_neural_network_bootstrap_test_performance\n",
    "\n",
    "\n",
    "\n",
    "        bootstrapped_stats.append(\n",
    "        {\n",
    "            'Sample ROC': ROC_AUC_neural_network_bootstrap_sample_performance,\n",
    "            'Test ROC': ROC_AUC_neural_network_bootstrap_test_performance,\n",
    "            'Optimisation': optimism\n",
    "        }\n",
    "       )\n",
    "\n",
    "\n",
    "bootstrapped_stats = pd.DataFrame(bootstrapped_stats)\n",
    "print(bootstrapped_stats.head())\n",
    "## Step 3: Get average optimization\n",
    "\n",
    "average_optimisation = bootstrapped_stats[\"Optimisation\"].mean() \n",
    "\n",
    "## Step 4: Get optimization-corrected performance\n",
    "\n",
    "optimization_corrected_performance = ROC_AUC_neural_network_apparent - average_optimisation ##\n",
    "\n",
    "print(optimization_corrected_performance)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "## Step 3: Get average optimization\n",
    "\n",
    "average_optimisation_nn = bootstrapped_stats[\"Optimisation\"].mean() \n",
    "\n",
    "## Step 4: Get optimization-corrected performance\n",
    "\n",
    "optimization_corrected_performance_nn = ROC_AUC_neural_network - average_optimisation_nn ##\n",
    "\n",
    "print(optimization_corrected_performance_nn)\n",
    "\n",
    "## get CI \n",
    "\n",
    "Bootstrap_CI_nn = (1 - 0.25)*bootstrapped_stats[\"Optimisation\"].quantile(q = 1)\n",
    "Upper_bootstrap_CI = optimization_corrected_performance_nn + Bootstrap_CI_nn\n",
    "Lower_bootstrap_CI = optimization_corrected_performance_nn - Bootstrap_CI_nn\n",
    "\n",
    "print(Upper_bootstrap_CI, Lower_bootstrap_CI)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "RandomOverSampler(sampling_strategy=0.5)\n",
      "0.7039313312560352\n"
     ]
    }
   ],
   "source": [
    "## Try oversampling (randomly duplicates examples in the minority class). Very imbalanced dataset, as CIPRO + is a strong minority \n",
    "## USE NO DROP DATA W/ DELTA REGION\n",
    "from imblearn.over_sampling import RandomOverSampler\n",
    "oversample = RandomOverSampler(sampling_strategy = 0.5)\n",
    "X = CIP_data[['MSM','MSMW', 'MSW', 'Oth/Unk/Missing','Midwest','Northeast', 'Southeast', 'Southwest', 'West', 'PREV_CLINIC', 'PREV_REGION', 'DELTA_REGION']]\n",
    "y = CIP_data['Susceptible']\n",
    "X_over, y_over = oversample.fit_resample(X,y)\n",
    "print(oversample)\n",
    "neural_network_model = MLPClassifier(solver = 'lbfgs', activation = 'tanh', max_iter = 1000 ,hidden_layer_sizes= 4, random_state=10)\n",
    "\n",
    "neural_network_fit = neural_network_model.fit(X_over,y_over)\n",
    "y_predict_nn = neural_network_fit.predict(X)\n",
    "\n",
    "ROC_AUC_neural_network = metrics.roc_auc_score(y, y_predict_nn)\n",
    "\n",
    "print(ROC_AUC_neural_network)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/model_selection/_validation.py:378: FitFailedWarning: \n",
      "200 fits failed out of a total of 1000.\n",
      "The score on these train-test partitions for these parameters will be set to nan.\n",
      "If these failures are not expected, you can try to debug them by setting error_score='raise'.\n",
      "\n",
      "Below are more details about the failures:\n",
      "--------------------------------------------------------------------------------\n",
      "149 fits failed with the following error:\n",
      "Traceback (most recent call last):\n",
      "  File \"/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/model_selection/_validation.py\", line 686, in _fit_and_score\n",
      "    estimator.fit(X_train, y_train, **fit_params)\n",
      "  File \"/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py\", line 740, in fit\n",
      "    self._validate_params()\n",
      "  File \"/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/base.py\", line 570, in _validate_params\n",
      "    validate_parameter_constraints(\n",
      "  File \"/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/utils/_param_validation.py\", line 97, in validate_parameter_constraints\n",
      "    raise InvalidParameterError(\n",
      "sklearn.utils._param_validation.InvalidParameterError: The 'solver' parameter of MLPClassifier must be a str among {'sgd', 'lbfgs', 'adam'}. Got 'sdg' instead.\n",
      "\n",
      "--------------------------------------------------------------------------------\n",
      "29 fits failed with the following error:\n",
      "Traceback (most recent call last):\n",
      "  File \"/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/model_selection/_validation.py\", line 686, in _fit_and_score\n",
      "    estimator.fit(X_train, y_train, **fit_params)\n",
      "  File \"/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py\", line 740, in fit\n",
      "    self._validate_params()\n",
      "  File \"/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/base.py\", line 570, in _validate_params\n",
      "    validate_parameter_constraints(\n",
      "  File \"/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/utils/_param_validation.py\", line 97, in validate_parameter_constraints\n",
      "    raise InvalidParameterError(\n",
      "sklearn.utils._param_validation.InvalidParameterError: The 'solver' parameter of MLPClassifier must be a str among {'lbfgs', 'adam', 'sgd'}. Got 'sdg' instead.\n",
      "\n",
      "--------------------------------------------------------------------------------\n",
      "18 fits failed with the following error:\n",
      "Traceback (most recent call last):\n",
      "  File \"/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/model_selection/_validation.py\", line 686, in _fit_and_score\n",
      "    estimator.fit(X_train, y_train, **fit_params)\n",
      "  File \"/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py\", line 740, in fit\n",
      "    self._validate_params()\n",
      "  File \"/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/base.py\", line 570, in _validate_params\n",
      "    validate_parameter_constraints(\n",
      "  File \"/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/utils/_param_validation.py\", line 97, in validate_parameter_constraints\n",
      "    raise InvalidParameterError(\n",
      "sklearn.utils._param_validation.InvalidParameterError: The 'solver' parameter of MLPClassifier must be a str among {'adam', 'lbfgs', 'sgd'}. Got 'sdg' instead.\n",
      "\n",
      "--------------------------------------------------------------------------------\n",
      "4 fits failed with the following error:\n",
      "Traceback (most recent call last):\n",
      "  File \"/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/model_selection/_validation.py\", line 686, in _fit_and_score\n",
      "    estimator.fit(X_train, y_train, **fit_params)\n",
      "  File \"/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py\", line 740, in fit\n",
      "    self._validate_params()\n",
      "  File \"/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/base.py\", line 570, in _validate_params\n",
      "    validate_parameter_constraints(\n",
      "  File \"/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/utils/_param_validation.py\", line 97, in validate_parameter_constraints\n",
      "    raise InvalidParameterError(\n",
      "sklearn.utils._param_validation.InvalidParameterError: The 'solver' parameter of MLPClassifier must be a str among {'lbfgs', 'sgd', 'adam'}. Got 'sdg' instead.\n",
      "\n",
      "  warnings.warn(some_fits_failed_message, FitFailedWarning)\n",
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/model_selection/_search.py:953: UserWarning: One or more of the test scores are non-finite: [0.80205189 0.80343402 0.80206016 0.8064685  0.77342634        nan\n",
      " 0.8064685         nan 0.80416979 0.80206016]\n",
      "  warnings.warn(\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Best Score: 0.8064685015234078\n",
      "Best Hyperparameters: {'solver': 'lbfgs', 'learning_rate': 'constant', 'hidden_layer_sizes': 12, 'alpha': 1.291549665014884, 'activation': 'relu'}\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Users/rem76/miniconda3/envs/GISP_init/lib/python3.10/site-packages/sklearn/neural_network/_multilayer_perceptron.py:536: ConvergenceWarning: lbfgs failed to converge (status=1):\n",
      "STOP: TOTAL NO. of ITERATIONS REACHED LIMIT.\n",
      "\n",
      "Increase the number of iterations (max_iter) or scale the data as shown in:\n",
      "    https://scikit-learn.org/stable/modules/preprocessing.html\n",
      "  self.n_iter_ = _check_optimize_result(\"lbfgs\", opt_res, self.max_iter)\n"
     ]
    }
   ],
   "source": [
    "### Try hyperparameter tuning with random oversampling \n",
    "from sklearn.model_selection import RepeatedStratifiedKFold\n",
    "from scipy.stats import loguniform\n",
    "from sklearn.model_selection import RandomizedSearchCV\n",
    "space = dict()\n",
    "space['activation'] = ['tanh', 'relu']\n",
    "space['solver'] = ['sdg', 'adam', 'lbfgs']\n",
    "space['alpha'] = np.logspace(-1, 1, 10)\n",
    "space['learning_rate'] = ['constant','adaptive']\n",
    "space['hidden_layer_sizes'] = [(4), (6),(8), (12)]\n",
    "\n",
    "cv = RepeatedStratifiedKFold(n_splits=10, n_repeats=10, random_state=1)\n",
    "\n",
    "\n",
    "model = MLPClassifier(random_state=10, max_iter = 3000)\n",
    "\n",
    "X = CIP_data[['MSM','MSMW', 'MSW', 'Oth/Unk/Missing','Midwest','Northeast', 'Southeast', 'Southwest', 'West', 'PREV_CLINIC', 'PREV_REGION', 'DELTA_REGION']]\n",
    "y = CIP_data['Susceptible']\n",
    "oversample = RandomOverSampler(sampling_strategy = 0.5)\n",
    "X, y = oversample.fit_resample(X,y)\n",
    "model_fit = model.fit(X, y)\n",
    "\n",
    "search = RandomizedSearchCV(model, space, scoring='roc_auc', n_jobs=-1, cv=cv, random_state=1)\n",
    "result = search.fit(X, y)\n",
    "# summarize result\n",
    "print('Best Score: %s' % result.best_score_)\n",
    "print('Best Hyperparameters: %s' % result.best_params_)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Try bootstrapping with oversampling and hyperparameters\n",
    "#Best Score: 0.8019328927793805\n",
    "#Best Hyperparameters: {'solver': 'lbfgs', 'learning_rate': 'constant', 'hidden_layer_sizes': 12, 'alpha': 1.291549665014884, 'activation': 'relu'}\n",
    "#bootstrap data\n",
    "n_iterations = 100\n",
    "bootstrapped_stats = []\n",
    "\n",
    "#1. Create model using all data and get ROC_AUC (\"ROC_AUC_neural_network\")\n",
    "X = CIP_data[['MSMW', 'MSW', 'Oth/Unk/Missing','Northeast', 'Southeast', 'Southwest', 'West', 'PREV_REGION', 'PREV_CLINIC', 'DELTA_REGION']]\n",
    "y = CIP_data['Susceptible']\n",
    "X, y = oversample.fit_resample(X,y)\n",
    "model = MLPClassifier(solver = 'lbfgs', activation = 'tanh', max_iter = 3000 ,hidden_layer_sizes= 12, alpha = 1.291549665014884, random_state=10, learning_rate = 'adaptive' )\n",
    "#model = MLPClassifier(solver = 'lbfgs', activation = 'relu', max_iter = 1000 ,hidden_layer_sizes= 12, alpha = 1.291549665014884, random_state=10, learning_rate = 'constant' )\n",
    "\n",
    "model_fit = model.fit(X, y)\n",
    "y_predict = model_fit.predict(X)\n",
    "\n",
    "ROC_AUC_neural_network_apparent = metrics.roc_auc_score(y, y_predict)\n",
    "\n",
    "for i in range(n_iterations):\n",
    "       #2. (A) Sample all individuals w/replacement\n",
    "        sample = CIP_data.sample(frac = 1, replace=True) ##(a) sample n individuals with replacement\n",
    "        X_sample = sample[['MSMW', 'MSW', 'Oth/Unk/Missing','Northeast', 'Southeast', 'Southwest', 'West', 'PREV_REGION', 'PREV_CLINIC', 'DELTA_REGION']]\n",
    "        y_sample = sample['Susceptible']\n",
    "        X_sample, y_sample = oversample.fit_resample(X_sample,y_sample)\n",
    "\n",
    "\n",
    "       #  (B) Develop predictive model and find apparent performance\n",
    "        model_fit = model.fit(X_sample, y_sample)\n",
    "        y_sample_predict = model_fit.predict(X_sample)\n",
    "        ROC_AUC_neural_network_bootstrap_sample_performance = metrics.roc_auc_score(y_sample, y_sample_predict) \n",
    "\n",
    "       #  (C) Performance of predictive model on original sample (i.e. original population, X)\n",
    "        y_test_predict = model_fit.predict(X)\n",
    "        ROC_AUC_neural_network_bootstrap_test_performance = metrics.roc_auc_score(y, y_test_predict) ## 0.756384214489288\n",
    "      ### (D) Calculate optimisation by getting (B) - (D) \n",
    "        optimism = ROC_AUC_neural_network_bootstrap_sample_performance - ROC_AUC_neural_network_bootstrap_test_performance\n",
    "\n",
    "\n",
    "\n",
    "        bootstrapped_stats.append(\n",
    "        {\n",
    "            'Sample ROC': ROC_AUC_neural_network_bootstrap_sample_performance,\n",
    "            'Test ROC': ROC_AUC_neural_network_bootstrap_test_performance,\n",
    "            'Optimisation': optimism\n",
    "        }\n",
    "       )\n",
    "\n",
    "\n",
    "bootstrapped_stats = pd.DataFrame(bootstrapped_stats)\n",
    "print(bootstrapped_stats.head())\n",
    "## Step 3: Get average optimization\n",
    "\n",
    "average_optimisation = bootstrapped_stats[\"Optimisation\"].mean() \n",
    "\n",
    "## Step 4: Get optimization-corrected performance\n",
    "\n",
    "optimization_corrected_performance = ROC_AUC_neural_network_apparent - average_optimisation ##\n",
    "\n",
    "print(optimization_corrected_performance)\n",
    "\n",
    "\n",
    "Bootstrap_CI_nn = (1 - 0.025)*bootstrapped_stats[\"Optimisation\"].quantile(q = 1)\n",
    "Upper_bootstrap_CI = optimization_corrected_performance + Bootstrap_CI_nn\n",
    "Lower_bootstrap_CI = optimization_corrected_performance - Bootstrap_CI_nn\n",
    "\n",
    "\n",
    "print(Upper_bootstrap_CI, Lower_bootstrap_CI)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Try bootstrapping with oversampling and hyperparameters and threshold \n",
    "from sklearn.metrics import roc_curve, auc \n",
    "from sklearn.metrics import confusion_matrix\n",
    "\n",
    "#bootstrap data\n",
    "n_iterations = 10\n",
    "bootstrapped_stats = []\n",
    "\n",
    "#1. Create model using all data and get ROC_AUC (\"ROC_AUC_neural_network\")\n",
    "X = CIP_data[['MSMW', 'MSW', 'Oth/Unk/Missing','Northeast', 'Southeast', 'Southwest', 'West', 'PREV_REGION', 'PREV_CLINIC', 'DELTA_REGION']]\n",
    "y = CIP_data['Susceptible']\n",
    "X, y = oversample.fit_resample(X,y)\n",
    "model = MLPClassifier(solver = 'lbfgs', activation = 'tanh', max_iter = 1000 ,hidden_layer_sizes= 12, random_state=10, learning_rate = 'adaptive' )\n",
    "model = MLPClassifier(solver = 'lbfgs', activation = 'tanh', max_iter = 3000 ,hidden_layer_sizes= 12, alpha = 1.291549665014884, random_state=10, learning_rate = 'adaptive' )\n",
    "\n",
    "model_fit = model.fit(X, y)\n",
    "y_predict = model_fit.predict(X)\n",
    "\n",
    "ROC_AUC_neural_network_apparent = metrics.roc_auc_score(y, y_predict)\n",
    "\n",
    "##threshold\n",
    "threshold = 0.5\n",
    "for i in range(n_iterations):\n",
    "       #2. (A) Sample all individuals w/replacement\n",
    "        sample = CIP_data.sample(frac = 1, replace=True) ##(a) sample n individuals with replacement\n",
    "        X_sample = sample[['MSMW', 'MSW', 'Oth/Unk/Missing','Northeast', 'Southeast', 'Southwest', 'West', 'PREV_REGION', 'PREV_CLINIC', 'DELTA_REGION']]]\n",
    "        y_sample = sample['Susceptible']\n",
    "        X_sample, y_sample = oversample.fit_resample(X_sample,y_sample)\n",
    "\n",
    "\n",
    "       #  (B) Develop predictive model and find apparent performance\n",
    "        model_fit = model.fit(X_sample, y_sample)\n",
    "        y_sample_predict = model_fit.predict(X_sample)\n",
    "        y_sample_predict_prob = model_fit.predict_proba(X_sample)\n",
    "\n",
    "        ROC_AUC_neural_network_bootstrap_sample_performance = metrics.roc_auc_score(y_sample, y_sample_predict) \n",
    "\n",
    "       #  (C) Performance of predictive model on original sample (i.e. original population, X)\n",
    "        y_test_predict = model_fit.predict(X)\n",
    "        y_test_predict_prob = model_fit.predict_proba(X)\n",
    "\n",
    "        ROC_AUC_neural_network_bootstrap_test_performance = metrics.roc_auc_score(y, y_test_predict) ## 0.756384214489288\n",
    "      ### (D) Calculate optimisation by getting (B) - (D) \n",
    "        optimism = ROC_AUC_neural_network_bootstrap_sample_performance - ROC_AUC_neural_network_bootstrap_test_performance\n",
    "     \n",
    "      ### (i) Calculate sensitivity and specificity \n",
    "        if threshold is not None:\n",
    "           y_sample_predict = np.where(y_sample_predict_prob[:, 1] > threshold, 1, 0)\n",
    "           y_test_predict = np.where(y_test_predict_prob[:, 1] > threshold, 1, 0)\n",
    "\n",
    "        tn_test, fp_test, fn_test, tp_test = confusion_matrix(y_true=y, y_pred=y_test_predict).ravel()\n",
    "        tn_sample, fp_sample, fn_sample, tp_sample = confusion_matrix(y_true=y_sample, y_pred=y_sample_predict).ravel()\n",
    "\n",
    "        sensitivity_test = tp_test / (tp_test + fn_test)\n",
    "        specificity_test = tn_test / (tn_test + fp_test)\n",
    "\n",
    "        sensitivity_sample = tp_sample/ (tp_sample + fn_sample)\n",
    "        specificity_sample = tn_sample / (tn_sample + fp_sample)\n",
    "\n",
    "        fpr_test, tpr_test, threshold_test = roc_curve(y, y_test_predict_prob[:, 1], drop_intermediate=False)\n",
    "        fpr_sample, tpr_sample, threshold_sample = roc_curve(y_sample, y_sample_predict_prob[:, 1], drop_intermediate=False)\n",
    "\n",
    "\n",
    "        #ROC_AUC_nn_test = roc_auc = auc(fpr_test, tpr_test)\n",
    "        #ROC_AUC_nn_sample = roc_auc = auc(fpr_sample, tpr_sample)\n",
    "\n",
    "\n",
    "        bootstrapped_stats.append(\n",
    "        {\n",
    "            'Sample ROC': ROC_AUC_neural_network_bootstrap_sample_performance,\n",
    "            'Test ROC': ROC_AUC_neural_network_bootstrap_test_performance,\n",
    "            'Optimisation': optimism,\n",
    "            'Sample Sensitivity': sensitivity_sample,\n",
    "            'Sample Specificity': specificity_sample,\n",
    "            'Test Sensitivity': sensitivity_test,\n",
    "            'Test Specificity': specificity_test,\n",
    "            'Sample FPR': fpr_sample,\n",
    "            'Sample TPR': tpr_sample,\n",
    "            'Test FPR': fpr_test,\n",
    "            'Test TPR': tpr_test\n",
    "        }\n",
    "       )\n",
    "\n",
    "\n",
    "bootstrapped_stats = pd.DataFrame(bootstrapped_stats)\n",
    "print(bootstrapped_stats.head())\n",
    "## Step 3: Get average optimization\n",
    "\n",
    "average_optimisation = bootstrapped_stats[\"Optimisation\"].mean() \n",
    "\n",
    "## Step 4: Get optimization-corrected performance\n",
    "\n",
    "optimization_corrected_performance = ROC_AUC_neural_network_apparent - average_optimisation ##\n",
    "\n",
    "print(optimization_corrected_performance)\n",
    "\n",
    "\n",
    "Bootstrap_CI_nn = (1 - 0.025)*bootstrapped_stats[\"Optimisation\"].quantile(q = 1)\n",
    "Upper_bootstrap_CI = optimization_corrected_performance + Bootstrap_CI_nn\n",
    "Lower_bootstrap_CI = optimization_corrected_performance - Bootstrap_CI_nn\n",
    "\n",
    "\n",
    "print(Upper_bootstrap_CI, Lower_bootstrap_CI)\n",
    "\n",
    "#### sensitivity "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#### Response of sensitivity and specificity to threshold\n",
    "\n",
    "from sklearn.metrics import roc_curve, auc \n",
    "from sklearn.metrics import confusion_matrix\n",
    "\n",
    "#bootstrap data\n",
    "n_iterations = 100\n",
    "bootstrapped_stats = []\n",
    "\n",
    "#1. Create model using all data and get ROC_AUC (\"ROC_AUC_neural_network\")\n",
    "\n",
    "model = MLPClassifier(solver = 'lbfgs', activation = 'tanh', max_iter = 2000 ,hidden_layer_sizes= 12, random_state=10, learning_rate = 'adaptive' )\n",
    "model = MLPClassifier(solver = 'lbfgs', activation = 'tanh', max_iter = 3000 ,hidden_layer_sizes= 12, alpha = 1.291549665014884, random_state=10, learning_rate = 'adaptive' )\n",
    "\n",
    "\n",
    "##threshold\n",
    "threshold_seq = [0,0.05,0.1,0.15,0.2,0.25,0.3,0.35,0.4,0.45,0.5,0.55,0.6,0.65,0.7,0.75,0.8,0.85,0.9,0.95,1]\n",
    "sensitivity_sample_table = []\n",
    "specificity_sample_table = []\n",
    "sensitivity_test_table = []\n",
    "specificity_test_table = []\n",
    "sensitivity_test_table_975 = []\n",
    "sensitivity_test_table_025 = []\n",
    "specificity_test_table_975 = []\n",
    "specificity_test_table_025 = []\n",
    "for threshold in threshold_seq:\n",
    "  print(threshold)\n",
    "  bootstrapped_stats = []\n",
    "  #1. Create model using all data and get ROC_AUC (\"ROC_AUC_neural_network\") - need to do this for each threshold\n",
    "\n",
    "  X = CIP_data[['MSMW', 'MSW', 'Oth/Unk/Missing','Northeast', 'Southeast', 'Southwest', 'West', 'PREV_REGION', 'PREV_CLINIC', 'DELTA_REGION']]\n",
    "  y = CIP_data['Susceptible']\n",
    "  X, y = oversample.fit_resample(X,y)\n",
    "  model_fit = model.fit(X, y)\n",
    "  y_predict = model_fit.predict(X)\n",
    "  y_predict_proba = model_fit.predict_proba(X)\n",
    " \n",
    "  y_predict = np.where(y_predict_proba[:, 1] > threshold, 1, 0)\n",
    "  tn_apparent , fp_apparent, fn_apparent, tp_apparent = confusion_matrix(y_true=y, y_pred=y_predict).ravel()\n",
    "\n",
    "  sensitivity_apparent = tp_apparent / (tp_apparent  + fn_apparent )\n",
    "  specificity_apparent  = tn_apparent / (tn_apparent + fp_apparent )\n",
    "  for i in range(n_iterations):\n",
    "        #2. (A) Sample all individuals w/replacement\n",
    "          sample = CIP_data.sample(frac = 1, replace=True) ##(a) sample n individuals with replacement\n",
    "          X_sample = sample[['MSMW', 'MSW', 'Oth/Unk/Missing','Northeast', 'Southeast', 'Southwest', 'West', 'PREV_REGION', 'PREV_CLINIC', 'DELTA_REGION']]\n",
    "          y_sample = sample['Susceptible']\n",
    "          X_sample, y_sample = oversample.fit_resample(X_sample,y_sample)\n",
    "\n",
    "\n",
    "        #  (B) Develop predictive model and find apparent performance\n",
    "          model_fit = model.fit(X_sample, y_sample)\n",
    "          y_sample_predict = model_fit.predict(X_sample)\n",
    "          y_sample_predict_prob = model_fit.predict_proba(X_sample)\n",
    "\n",
    "          ROC_AUC_neural_network_bootstrap_sample_performance = metrics.roc_auc_score(y_sample, y_sample_predict) \n",
    "\n",
    "        #  (C) Performance of predictive model on original sample (i.e. original population, X)\n",
    "          y_test_predict = model_fit.predict(X)\n",
    "          y_test_predict_prob = model_fit.predict_proba(X)\n",
    "\n",
    "          ROC_AUC_neural_network_bootstrap_test_performance = metrics.roc_auc_score(y, y_test_predict) ## 0.756384214489288\n",
    "        ### (D) Calculate optimisation by getting (B) - (D) \n",
    "          optimism = ROC_AUC_neural_network_bootstrap_sample_performance - ROC_AUC_neural_network_bootstrap_test_performance\n",
    "      \n",
    "        ### (i) Calculate sensitivity and specificity \n",
    "          \n",
    "          y_sample_predict = np.where(y_sample_predict_prob[:, 1] > threshold, 1, 0)\n",
    "          y_test_predict = np.where(y_test_predict_prob[:, 1] > threshold, 1, 0)\n",
    "\n",
    "          tn_test, fp_test, fn_test, tp_test = confusion_matrix(y_true=y, y_pred=y_test_predict).ravel()\n",
    "          tn_sample, fp_sample, fn_sample, tp_sample = confusion_matrix(y_true=y_sample, y_pred=y_sample_predict).ravel()\n",
    "\n",
    "          sensitivity_test = tp_test / (tp_test + fn_test)\n",
    "          specificity_test = tn_test / (tn_test + fp_test)\n",
    "\n",
    "          sensitivity_sample = tp_sample/ (tp_sample + fn_sample)\n",
    "          specificity_sample = tn_sample / (tn_sample + fp_sample)\n",
    "\n",
    "          fpr_test, tpr_test, threshold_test = roc_curve(y, y_test_predict_prob[:, 1], drop_intermediate=False)\n",
    "          fpr_sample, tpr_sample, threshold_sample = roc_curve(y_sample, y_sample_predict_prob[:, 1], drop_intermediate=False)\n",
    "          \n",
    "          sensitivity_optimism = sensitivity_sample - sensitivity_test\n",
    "          specificity_optimism = specificity_sample - specificity_test\n",
    "\n",
    "          bootstrapped_stats.append(\n",
    "        {\n",
    "            'Sample ROC': ROC_AUC_neural_network_bootstrap_sample_performance,\n",
    "            'Test ROC': ROC_AUC_neural_network_bootstrap_test_performance,\n",
    "            'Optimisation': optimism,\n",
    "            'Sample Sensitivity': sensitivity_sample,\n",
    "            'Sample Specificity': specificity_sample,\n",
    "            'Test Sensitivity': sensitivity_test,\n",
    "            'Test Specificity': specificity_test,\n",
    "            'Sensitivity Optimisation': sensitivity_optimism,\n",
    "            'Specificity Optimisation': specificity_optimism,\n",
    "            'Sample FPR': fpr_sample,\n",
    "            'Sample TPR': tpr_sample,\n",
    "            'Test FPR': fpr_test,\n",
    "            'Test TPR': tpr_test\n",
    "        })\n",
    "  bootstrapped_stats = pd.DataFrame(bootstrapped_stats)\n",
    "  specificity_sample_optimism = bootstrapped_stats[\"Specificity Optimisation\"].mean() \n",
    "  sensitivity_sample_optimism = bootstrapped_stats[\"Sensitivity Optimisation\"].mean() \n",
    "\n",
    "  specificity_sample_optimised = specificity_apparent - specificity_sample_optimism ##\n",
    "  sensitivity_sample_optimised = sensitivity_apparent - sensitivity_sample_optimism ##\n",
    "\n",
    "  specificity_sample_95CI_lower = specificity_sample_optimised - bootstrapped_stats[\"Specificity Optimisation\"].quantile(0.025)\n",
    "  specificity_sample_95CI_upper = specificity_sample_optimised + bootstrapped_stats[\"Specificity Optimisation\"].quantile(0.975)\n",
    "\n",
    "  sensitivity_sample_95CI_lower= sensitivity_sample_optimised - bootstrapped_stats[\"Sensitivity Optimisation\"].quantile(0.025)\n",
    "  sensitivity_sample_95CI_upper = sensitivity_sample_optimised + bootstrapped_stats[\"Sensitivity Optimisation\"].quantile(0.975)\n",
    "\n",
    "  #sensitivity_sample_table.append(sensitivity_sample_optimised)\n",
    "  #specificity_sample_table.append(specificity_sample_optimised)\n",
    "  sensitivity_test_table.append(sensitivity_sample_optimised)\n",
    "  specificity_test_table.append(specificity_sample_optimised)\n",
    "  sensitivity_test_table_975.append(sensitivity_sample_95CI_upper)\n",
    "  sensitivity_test_table_025.append(sensitivity_sample_95CI_lower)\n",
    "  specificity_test_table_975.append(specificity_sample_95CI_upper)\n",
    "  specificity_test_table_025.append(specificity_sample_95CI_lower)\n",
    "\n",
    "          \n",
    "        \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#print(specificity_test_table)\n",
    "#print(sensitivity_test_table)\n",
    "print(specificity_test_table)\n",
    "\n",
    "print(specificity_test_table_975)\n",
    "\n",
    "#specificity_test_table = pd.DataFrame(specificity_test_table)\n",
    "#sensitivity_test_table = pd.DataFrame(sensitivity_test_table)\n",
    "#specificity_test_table.haed()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sensitivity_test_df = pd.DataFrame(sensitivity_test_table)\n",
    "sensitivity_test_df_975 = pd.DataFrame(sensitivity_test_table_975)\n",
    "sensitivity_test_df_025 = pd.DataFrame(sensitivity_test_table_025)\n",
    "\n",
    "specificity_test_df_975 = pd.DataFrame(specificity_test_table_975)\n",
    "specificity_test_df_025 = pd.DataFrame(specificity_test_table_025)\n",
    "\n",
    "threshold_seq = [0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1]\n",
    "plt.plot(threshold_seq, sensitivity_test_table, color = \"black\", label = \"Sensitivity\")\n",
    "plt.plot(threshold_seq, specificity_test_table, color = \"blue\", label = \"Specificity\")\n",
    "#error\n",
    "plt.plot(threshold_seq, sensitivity_test_df_975, color = \"slategrey\", label = \"\", linewidth=0.5)\n",
    "plt.plot(threshold_seq, sensitivity_test_df_025, color = \"slategrey\", label = \"\", linewidth=0.5)\n",
    "plt.plot(threshold_seq, specificity_test_df_975, color = \"slategrey\", label = \"\", linewidth=0.5)\n",
    "plt.plot(threshold_seq, specificity_test_df_025, color = \"slategrey\", label = \"\", linewidth=0.5)\n",
    "\n",
    "\n",
    "\n",
    "#plt.errorbar(threshold_seq, specificity_test_table, yerr = y_err_sensitivity, ecolor=\"lightblue\")\n",
    "\n",
    "plt.xlabel(\"Threshold\", fontsize=12)\n",
    "plt.ylabel(\"Performance metric\", fontsize=12)\n",
    "plt.xticks(threshold_seq, fontsize=12)\n",
    "plt.rcParams['figure.dpi'] = 300\n",
    "plt.legend()"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "GISP_init",
   "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.10.9 (main, Jan 11 2023, 09:18:18) [Clang 14.0.6 ]"
  },
  "orig_nbformat": 4,
  "vscode": {
   "interpreter": {
    "hash": "37b4cc77837642ee25c4eb6578aebe03e17eb3bb59efdde49edbbc888dbc418f"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}