{"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":"#  Implementing Logistic Regression from scratch\n","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2022-10-16T05:29:35.619607Z","iopub.execute_input":"2022-10-16T05:29:35.620511Z","iopub.status.idle":"2022-10-16T05:29:35.635606Z","shell.execute_reply.started":"2022-10-16T05:29:35.620467Z","shell.execute_reply":"2022-10-16T05:29:35.634856Z"}}},{"cell_type":"code","source":"import numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\nfrom sklearn.model_selection import train_test_split\nimport math\nimport matplotlib.pyplot as plt\nfrom numpy import linalg as LA\n\n\n# Input data files are available in the read-only \"../input/\" directory\n# For example, running this (by clicking run or pressing Shift+Enter) will list all files under the input directory\n\nimport os\nfor dirname, _, filenames in os.walk('/kaggle/input'):\n    for filename in filenames:\n        print(os.path.join(dirname, filename))\n\n\n\n# You can write up to 20GB to the current directory (/kaggle/working/) that gets preserved as output when you create a version using \"Save & Run All\" \n# You can also write temporary files to /kaggle/temp/, but they won't be saved outside of the current session\n","metadata":{"execution":{"iopub.status.busy":"2022-12-20T07:23:40.313555Z","iopub.execute_input":"2022-12-20T07:23:40.313982Z","iopub.status.idle":"2022-12-20T07:23:40.326137Z","shell.execute_reply.started":"2022-12-20T07:23:40.313950Z","shell.execute_reply":"2022-12-20T07:23:40.324569Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"data = pd.read_csv(\"/kaggle/input/heart-disease-prediction-using-logistic-regression/framingham.csv\")","metadata":{"execution":{"iopub.status.busy":"2022-12-20T07:23:40.328619Z","iopub.execute_input":"2022-12-20T07:23:40.329000Z","iopub.status.idle":"2022-12-20T07:23:40.361906Z","shell.execute_reply.started":"2022-12-20T07:23:40.328968Z","shell.execute_reply":"2022-12-20T07:23:40.360124Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"data.head()","metadata":{"execution":{"iopub.status.busy":"2022-12-20T07:23:40.363764Z","iopub.execute_input":"2022-12-20T07:23:40.364609Z","iopub.status.idle":"2022-12-20T07:23:40.387722Z","shell.execute_reply.started":"2022-12-20T07:23:40.364568Z","shell.execute_reply":"2022-12-20T07:23:40.386295Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"data.shape\ndata.isnull().sum()","metadata":{"execution":{"iopub.status.busy":"2022-12-20T07:23:40.389945Z","iopub.execute_input":"2022-12-20T07:23:40.390346Z","iopub.status.idle":"2022-12-20T07:23:40.400944Z","shell.execute_reply.started":"2022-12-20T07:23:40.390309Z","shell.execute_reply":"2022-12-20T07:23:40.399945Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Also I *very naively*, managed to fill the NA values in the dataset!","metadata":{}},{"cell_type":"code","source":"mean_value_glucose = data['glucose'].mean()\ndata['glucose'].fillna(value=mean_value_glucose, inplace=True)\ndata['cigsPerDay'].fillna(value=0.0, inplace=True)\nnp.random.seed(3)\ndata['education'].fillna(value=np.random.randint(1,5), inplace=True)\nmean_value_chol = data['totChol'].mean()\ndata['totChol'].fillna(value=mean_value_chol, inplace=True)\ndata['BPMeds'].fillna(value=np.random.randint(0,2), inplace=True)\ndata['BMI'].fillna(value=np.random.randint(0,2), inplace=True)\nmean_value_rate = data['heartRate'].mean()\ndata['heartRate'].fillna(value=mean_value_rate, inplace=True)","metadata":{"execution":{"iopub.status.busy":"2022-12-20T07:23:40.402643Z","iopub.execute_input":"2022-12-20T07:23:40.403143Z","iopub.status.idle":"2022-12-20T07:23:40.421323Z","shell.execute_reply.started":"2022-12-20T07:23:40.403111Z","shell.execute_reply":"2022-12-20T07:23:40.420060Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#  Parsing out the target from training data.\n\ny = data['TenYearCHD']\nX = data.drop('TenYearCHD', axis = 1)\nX = X.to_numpy()\ny = y.to_numpy()\n","metadata":{"execution":{"iopub.status.busy":"2022-12-20T07:23:40.422973Z","iopub.execute_input":"2022-12-20T07:23:40.423466Z","iopub.status.idle":"2022-12-20T07:23:40.444229Z","shell.execute_reply.started":"2022-12-20T07:23:40.423410Z","shell.execute_reply":"2022-12-20T07:23:40.442058Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Since the data has its features ranging very differently from each other, we do a\nfavour to the gradient descent to help the training over the data set so that it \nconverges much faster and efficiently with a good learning rate.\n\nWe do feature scaling!\nWhy?\n\n1. The learning rate controls the size of the update to the parameters during training. This size of updation is same in this notebook for each of the paramaters of the features.\n\n2. What varies is the X_j term(the jth feature value of ith training set). The set of features with a higher range need the training of their parameters such that they have a low value as a small change in the parameter of such a feature would have a great impact on the estimated prediction and hence also the cost function! \n\n3. After learning the parameters from the model, when we predict the values, from the features we have not seen before. Given a new X value, we must first scale x using the min_vec and max_vec that we had previously computed from the training set.\n\n\nHow?\n\n1. When we scale features, the scaling constants like min and max value, standard deviation etc, come from only the training set.\n\n\n","metadata":{}},{"cell_type":"code","source":"def scale_features(X):\n    \"\"\"\n    \n    Returns:\n      X_scale (ndarray (m,n)): input scaled by column\n      min_ (ndarray (n,))     : minimum of each feature\n      max_ (ndarray (n,))  : maximum of each feature\n      \n    \"\"\"\n    min_vec = X.min(axis = 0)\n    max_vec = X.max(axis = 0)\n    X = (X - min_vec) / (max_vec - min_vec)\n    return X, min_vec, max_vec","metadata":{"execution":{"iopub.status.busy":"2022-12-20T07:23:40.446239Z","iopub.execute_input":"2022-12-20T07:23:40.446728Z","iopub.status.idle":"2022-12-20T07:23:40.459127Z","shell.execute_reply.started":"2022-12-20T07:23:40.446688Z","shell.execute_reply":"2022-12-20T07:23:40.457592Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.25, random_state=1)\nprint(\"X_train.shape\", X_train.shape, \"y_train.shape\", y_train.shape)\nprint(\"X_test.shape\", X_test.shape, \"y_test.shape\", y_test.shape)\n","metadata":{"execution":{"iopub.status.busy":"2022-12-20T07:23:40.461048Z","iopub.execute_input":"2022-12-20T07:23:40.462319Z","iopub.status.idle":"2022-12-20T07:23:40.479557Z","shell.execute_reply.started":"2022-12-20T07:23:40.462257Z","shell.execute_reply":"2022-12-20T07:23:40.477985Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"X_train, min_vec, max_vec = scale_features(X_train)\nprint(X_train)","metadata":{"execution":{"iopub.status.busy":"2022-12-20T07:23:40.481037Z","iopub.execute_input":"2022-12-20T07:23:40.482004Z","iopub.status.idle":"2022-12-20T07:23:40.491544Z","shell.execute_reply.started":"2022-12-20T07:23:40.481957Z","shell.execute_reply":"2022-12-20T07:23:40.490560Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def sigmoid(z) :\n    return 1 / (1 + np.exp(-z))\n","metadata":{"execution":{"iopub.status.busy":"2022-12-20T07:23:40.566598Z","iopub.execute_input":"2022-12-20T07:23:40.567878Z","iopub.status.idle":"2022-12-20T07:23:40.573467Z","shell.execute_reply.started":"2022-12-20T07:23:40.567832Z","shell.execute_reply":"2022-12-20T07:23:40.571958Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# COST FUNCTION FOR LOGISTIC REGRESSION :\n\n\ndef compute_cost(X, y, w, b) :\n    m, n = X.shape\n    \n    # loop over each example\n    loss = 0\n    for example in range(m):\n        z = 0\n        for feature in range(n):\n            z += X[example, feature] * w[feature]\n        z += b\n        f_wb = sigmoid(z)\n        loss += - y[example] * np.log(f_wb) - (1 - y[example]) * np.log(1 - f_wb)\n    return (loss / m)\n\n\nnp.random.seed(1)\ncompute_cost(X_train, y_train, 0.01 * np.random.rand(15), 1)\n","metadata":{"execution":{"iopub.status.busy":"2022-12-20T07:23:40.575377Z","iopub.execute_input":"2022-12-20T07:23:40.576591Z","iopub.status.idle":"2022-12-20T07:23:40.640401Z","shell.execute_reply.started":"2022-12-20T07:23:40.576543Z","shell.execute_reply":"2022-12-20T07:23:40.638720Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"\nImplementing Gradient for Logistic Regression :","metadata":{}},{"cell_type":"code","source":"# # Vectorised implementation : \n\n# # def compute_gradient(X, y, w, b, lambda_ = None) :\n    \n# #     m, n = X.shape\n# #     dj_dw, dj_db = np.zeros(n), 0.\n    \n    \n# #     f_wb = sigmoid(np.dot(X, w) + b)           \n# #     err = (f_wb - y).reshape(-1, 1)\n# #     print(err)\n# #     dj_dw = sum(err * X) / m\n# #     dj_db = sum(err) / m\n    \n# #     return dj_dw, dj_db\n\n# # Core Implementation :\n\ndef compute_gradient(X, y, w, b, lambda_ = None) :\n    m, n = X.shape\n    dj_dw, dj_db = np.zeros(n), 0.\n    \n    for eg in range(m) :   \n        \n        f_wb = sigmoid(np.dot(X[eg], w) + b)\n        err = f_wb - y[eg]\n        \n        for fr in range(n) :\n            dj_dw[fr] += err * X[eg, fr]\n        dj_db += err\n            \n    dj_dw /= m\n    dj_db /= m\n    return dj_db, dj_dw\n\ncompute_gradient(X_train, y_train, 0.01 * np.random.rand(15), 1, lambda_ = None)","metadata":{"execution":{"iopub.status.busy":"2022-12-20T07:23:40.642153Z","iopub.execute_input":"2022-12-20T07:23:40.642643Z","iopub.status.idle":"2022-12-20T07:23:40.729806Z","shell.execute_reply.started":"2022-12-20T07:23:40.642600Z","shell.execute_reply":"2022-12-20T07:23:40.728320Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"We would implement gradient descent. To verify that gradient descent is working properly, it would be great if we can see the value of the cost function decreasing with each step of descent.","metadata":{}},{"cell_type":"code","source":"\ndef gradient_descent(X, y, w_in, b_in, alpha, num_iters, cost_func, gradient_func, lambda_) :   \n         \n    for i in range(num_iters) :\n        \n        dj_db, dj_dw = gradient_func(X, y, w_in, b_in, lambda_)\n        w_in = w_in - alpha * dj_dw\n        b_in = b_in - alpha * dj_db\n        \n#         For 10 iterations equally spaced of total iters, print cost.\n        if i % math.ceil(num_iters / 10) == 0 or i == (num_iters - 1) :\n            print(f\"Iteration {i : 4}: Cost {cost_func(X, y, w_in, b_in) : 8.2f}\")\n            \n    return w_in, b_in   ","metadata":{"execution":{"iopub.status.busy":"2022-12-20T07:23:40.731604Z","iopub.execute_input":"2022-12-20T07:23:40.732492Z","iopub.status.idle":"2022-12-20T07:23:40.741499Z","shell.execute_reply.started":"2022-12-20T07:23:40.732418Z","shell.execute_reply":"2022-12-20T07:23:40.740161Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"* * **Running the gradient descent to learn the parameters for our dataset:**","metadata":{}},{"cell_type":"code","source":"\nnp.random.seed(1)\ninitial_w = 0.01* (np.random.rand(15) - 0.5)\nprint(initial_w)\ninitial_b = -2\n\n# Setting up gradient descent parameters.\n\n# Please try 10000 iteration to see how slow the descent is!\n# num_iters = 10000\n\n# Since gradient descent is too slow , trying with 1000 iterations\nnum_iters = 1000\n\n# Feel free to try with different learning rates:\nalpha = 0.5\n# alpha = 0.000001\n# Without Regularisation:\nlambda_ = 0\n\nw, b = gradient_descent(X_train, y_train, initial_w, initial_b, alpha, num_iters,compute_cost, compute_gradient, lambda_)","metadata":{"execution":{"iopub.status.busy":"2022-12-20T07:23:40.744527Z","iopub.execute_input":"2022-12-20T07:23:40.744928Z","iopub.status.idle":"2022-12-20T07:24:11.227282Z","shell.execute_reply.started":"2022-12-20T07:23:40.744896Z","shell.execute_reply":"2022-12-20T07:24:11.225890Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"The cost function decreases more efficiently than before which a larger value of alpha.\n\n\nIf above, the cost function was not decreasing as the gradient descents, then:\n\n1. implies that learning rate is too high relatively!\n2. Some data features need to be engineered?\n3. More data needs to be collected?\n\nLet us also compute traing error and test error on the learnt parameters:","metadata":{}},{"cell_type":"code","source":"compute_cost(X_train, y_train, w, b)","metadata":{"execution":{"iopub.status.busy":"2022-12-20T07:24:11.228593Z","iopub.execute_input":"2022-12-20T07:24:11.228937Z","iopub.status.idle":"2022-12-20T07:24:11.283965Z","shell.execute_reply.started":"2022-12-20T07:24:11.228904Z","shell.execute_reply":"2022-12-20T07:24:11.282490Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"We would now scale test data using normalised features!","metadata":{}},{"cell_type":"code","source":"X_test = (X_test - min_vec) / (max_vec - min_vec)\nprint(X_test)","metadata":{"execution":{"iopub.status.busy":"2022-12-20T07:24:11.286052Z","iopub.execute_input":"2022-12-20T07:24:11.286449Z","iopub.status.idle":"2022-12-20T07:24:11.293294Z","shell.execute_reply.started":"2022-12-20T07:24:11.286393Z","shell.execute_reply":"2022-12-20T07:24:11.291919Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"compute_cost(X_test, y_test, w, b)","metadata":{"execution":{"iopub.status.busy":"2022-12-20T07:24:11.294690Z","iopub.execute_input":"2022-12-20T07:24:11.295018Z","iopub.status.idle":"2022-12-20T07:24:11.328146Z","shell.execute_reply.started":"2022-12-20T07:24:11.294985Z","shell.execute_reply":"2022-12-20T07:24:11.327039Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Seems counter intuitive! The test error is lower than training error,and therefore\nthe essence of machine learning that is to predict for the unknown is disregarded by the model as the model does more good to predict than to train, it does even worse on the training set!\n\nWhat our model has learnt? If the test error is even lower than the train error, what was the point of training then?","metadata":{}},{"cell_type":"markdown","source":"Apart from above doubts, Lets predict the values by our model!","metadata":{}},{"cell_type":"code","source":"def predict(X, w, b) :\n    m, n = X.shape\n    p = np.zeros(m)    \n    for i in range(m) :\n        f_wb = sigmoid(np.dot(X[i], w) + b)\n        p[i] = 1 if f_wb > 0.38 else 0\n    return p","metadata":{"execution":{"iopub.status.busy":"2022-12-20T07:24:11.329027Z","iopub.execute_input":"2022-12-20T07:24:11.329330Z","iopub.status.idle":"2022-12-20T07:24:11.336765Z","shell.execute_reply.started":"2022-12-20T07:24:11.329301Z","shell.execute_reply":"2022-12-20T07:24:11.335010Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"p = predict(X_train, w,b)\nprint('Train Accuracy: %f'%(np.mean(p == y_train) * 100))\np = predict(X_test, w,b)\nprint('Test Accuracy: %f'%(np.mean(p == y_test) * 100))","metadata":{"execution":{"iopub.status.busy":"2022-12-20T07:24:11.338285Z","iopub.execute_input":"2022-12-20T07:24:11.338814Z","iopub.status.idle":"2022-12-20T07:24:11.363683Z","shell.execute_reply.started":"2022-12-20T07:24:11.338770Z","shell.execute_reply":"2022-12-20T07:24:11.362559Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}