{"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":"# Setup","metadata":{"_kg_hide-input":true}},{"cell_type":"code","source":"import pandas as pd\nimport numpy as np\n\nimport matplotlib.pyplot as plt\nplt.style.use(\"seaborn-whitegrid\")\nplt.rc(\"figure\", autolayout=True)\nplt.rc(\n    \"axes\",\n    labelweight=\"bold\",\n    labelsize=\"large\",\n    titleweight=\"bold\",\n    titlesize=14,\n    titlepad=10,\n)\n\nimport warnings\n\nwarnings.filterwarnings('ignore')\n\nimport xgboost as xgb\nimport sklearn\nfrom sklearn.metrics import accuracy_score,classification_report,\\\n                        confusion_matrix","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-08-01T11:32:40.311302Z","iopub.execute_input":"2022-08-01T11:32:40.311818Z","iopub.status.idle":"2022-08-01T11:32:41.423169Z","shell.execute_reply.started":"2022-08-01T11:32:40.311704Z","shell.execute_reply":"2022-08-01T11:32:41.422057Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"def plot_confusion_matrix(cm,\n                          target_names,\n                          title='Confusion matrix',\n                          cmap=None,\n                          normalize=True):\n    import matplotlib.pyplot as plt\n    import numpy as np\n    import itertools\n\n    accuracy = np.trace(cm) / float(np.sum(cm))\n    misclass = 1 - accuracy\n\n    if cmap is None:\n        cmap = plt.get_cmap('Blues')\n\n    plt.figure(figsize=(8, 6))\n    plt.imshow(cm, interpolation='nearest', cmap=cmap)\n    plt.title(title)\n    plt.colorbar()\n\n    if target_names is not None:\n        tick_marks = np.arange(len(target_names))\n        plt.xticks(tick_marks, target_names, rotation=45)\n        plt.yticks(tick_marks, target_names)\n\n    if normalize:\n        cm = cm.astype('float') / cm.sum(axis=1)[:, np.newaxis]\n\n\n    thresh = cm.max() / 1.5 if normalize else cm.max() / 2\n    for i, j in itertools.product(range(cm.shape[0]), range(cm.shape[1])):\n        if normalize:\n            plt.text(j, i, \"{:0.4f}\".format(cm[i, j]),\n                     horizontalalignment=\"center\",\n                     color=\"white\" if cm[i, j] > thresh else \"black\")\n        else:\n            plt.text(j, i, \"{:,}\".format(cm[i, j]),\n                     horizontalalignment=\"center\",\n                     color=\"white\" if cm[i, j] > thresh else \"black\")\n\n\n    plt.tight_layout()\n    plt.ylabel('True label')\n    plt.xlabel('Predicted label\\naccuracy={:0.4f}; misclass={:0.4f}'.format(accuracy, misclass))\n    plt.show()","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-08-01T11:32:41.964410Z","iopub.execute_input":"2022-08-01T11:32:41.964718Z","iopub.status.idle":"2022-08-01T11:32:41.982936Z","shell.execute_reply.started":"2022-08-01T11:32:41.964688Z","shell.execute_reply":"2022-08-01T11:32:41.981625Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Introduction\nXGBoost is a very powerful algorithm that has one of the best performance on tabular data for regression and classification tasks.\n\nIt has many important features:\n- easy to start;\n- good performance with default hyperparameters;\n- has own API, but a user can choose sklearn API to take all advantages of scikit-learn.\n- wide range of methods to prevent overfitting\n- visualisation of decision trees and SHAP values to get insights about features impact.\nAnd many more important things.\n\nXGBoost is designed to be an extensible library. And one of the important techniques is creating custom objective (loss) functions.\nExamples of a creating custom loss functions for regression tasks provided in XGBoost documentation and in various articles.\nSome references and helpful sources would be provided in the end of the notebook.\n\nThe topic of the notebook is a creation of a custom loss function and a respective metric for binary classification problems.","metadata":{}},{"cell_type":"markdown","source":"## Why to use custom loss function?\nWell, there may be situations when accuracy metric is not enough to get expected result. We may want to reduce false Negative (FN) or false positive (FP) rate. It may be very helpful when dataset is imbalanced and the result we seek belongs to a smaller class.\n\nExamples:\n- fraud detection. We want to reduce false negative samples (model treated transaction as usual when actually it is doubtful) to reduce risk that fraudulent transaction bypass unnoticed (even if false positive rate would grow up).\n- disease identification. It has the same logic we want to increase chances to catch the disease.\n\nOf course, in both cases we do not want to get extremely high FP rates, because it would increase resource usage to confirm that sample is not fraud/disease. So model optimization remains necessary.","metadata":{}},{"cell_type":"markdown","source":"# A little bit of math\nA standard objective function for a binary classification problem in XGBoost is \"objective\":\"binary:logistic\". The goal of the algortihm is to find weights of features that minimise logloss (binary crossentropy) function:\n\n$$ l(y,p)=-\\frac{1}{N}\\sum_{i=1}^{N} (y_{i}log(p_{i})+(1-y_{i})log(1-p_{i})) $$ [ammendment in response to Thomas Aarholt comment]\n\nwhere N - number of samples, y - true label, log - natural logarithm.\n\np_i - probability of the ith sample being 1 (sigmoid function):\n\n$$ p_{i}=\\frac{1}{1+e^{-z}}$$\n\nwhere m - number of features,\n\n$$ z=w_{0}+w_{1}x_{1}+...+w_{m}x_{m}$$\n\nz is linear regression where w - weights of features, x - feature value\n\nLet's plot logloss function to visualise its shape.","metadata":{}},{"cell_type":"code","source":"# plot logloss function (logloss for label 1, logloss for label 0, overall logloss)\n# plot logloss function:\n\n# create p and true labels distribution\np=np.linspace(1e-04,.9999,endpoint=True)\ny=np.random.binomial(n=1,p=.5,size=[50])\n\n# create three subplots\nplt.figure(figsize=(20,5))\nplt.subplot(131)\nplt.plot(p,-np.log(p))\nplt.plot(p,-np.log(p),lw=3)\nplt.xlabel('probability')\nplt.ylabel('loss: -log(p)')\nplt.title('Logloss for class 1')\n\nplt.subplot(132)\nplt.plot(p,-np.log(1-p))\nplt.plot(p,-np.log(1-p),lw=3)\nplt.xlabel('probability')\nplt.ylabel('loss: -log(1-p)')\nplt.title('Logloss for class 0')\n\nplt.subplot(133)\nplt.plot(p,-np.log(p)-np.log(1-p))\nplt.plot(p,-np.log(p)-np.log(1-p),lw=3)\nplt.xlabel('probability')\nplt.ylabel('loss: -log(p)-log(1-p)')\nplt.title('Cumulative logloss')","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-08-01T11:33:06.054740Z","iopub.execute_input":"2022-08-01T11:33:06.055716Z","iopub.status.idle":"2022-08-01T11:33:06.860577Z","shell.execute_reply.started":"2022-08-01T11:33:06.055659Z","shell.execute_reply":"2022-08-01T11:33:06.859546Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"\nThe first graph shows max error due to probability of predicted class being 1 is nearly 0% and actual class is 1. It is case of false negative.\n\nThe second graph is an opposite case. And when probability goes to 100%, it leads to false positive result.\n\nThe third graph is cumulative loss function. Logloss function is convex and has a single minimum.\n","metadata":{}},{"cell_type":"markdown","source":"## Calculation of first order and second order derivatives for general logloss function\nIt is needed to calculate first order gradient (gradient) and second order gradient (hessian) to update loss function, based on model predictions and observed data labels.\n\n### Important note\nVariable for differentiation is function z (linear regression). Thus chain rule should be applied:\n\n$$ \\frac{\\delta l}{\\delta z}=\\frac{\\delta l}{\\delta p}\\frac{\\delta p}{\\delta z} $$\n\nFinal elementwise derivatives are as follows:\n\n$$ gradient = p_{i} - y_{i} $$\n\n$$ hessian = p_{i}(1-p_{i}) $$\n\n\n## Introduction of custom logloss function and calculating its first order and second order derivatives\nLet's introduce a parameter beta in range (0,1) as follows: if beta>1 than model would be penalized for false positives, if beta<1 than false negatives.\n\nSo logloss function becomes:\n\n$$ l(y,p)=-\\frac{1}{N}\\sum_{i=1}^{N} (y_{i}log(p_{i})+\\beta(1-y_{i})log(1-p_{i})) $$ [ammendment in response to Thomas Aarholt comment]\n\n\n$$ gradient = p_{i}(\\beta - \\beta y_{i} + y_{i} ) - y_{i} $$\n$$ hessian = ( y_{i} + \\beta - \\beta y_{i} )p_{i}(1-p_{i}) $$\n\nAs you can see probability and targets are needed for inputs in user defined function.\n\nLet's plot custom loss functions with beta=0.4 and beta=2. Plots show the impact of beta to penalize FN (when beta=0.4) or FP (when beta=2)","metadata":{}},{"cell_type":"code","source":"plt.figure(figsize=(10,5))\nplt.subplot(121)\nplt.plot(p,-np.log(p)-0.4*np.log(1-p))\nplt.plot(p,-np.log(p)-0.4*np.log(1-p),lw=3)\nplt.xlabel('probability')\nplt.ylabel('loss: -log(p)-0.4*log(1-p)')\nplt.title('Cumulative logloss with beta 0.4')\n\nplt.subplot(122)\nplt.plot(p,-np.log(p)-2*np.log(1-p))\nplt.plot(p,-np.log(p)-2*np.log(1-p),lw=3)\nplt.xlabel('probability')\nplt.ylabel('loss: -log(p)-2*log(1-p)')\nplt.title('Cumulative logloss with beta 2')","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2022-08-01T11:34:42.246169Z","iopub.execute_input":"2022-08-01T11:34:42.246486Z","iopub.status.idle":"2022-08-01T11:34:42.771940Z","shell.execute_reply.started":"2022-08-01T11:34:42.246456Z","shell.execute_reply":"2022-08-01T11:34:42.770939Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Examples (finally)\nI have chosen Titanic dataset because it is approachable, imbalanced and easy to use.\nFirst, I have prepared a baseline model with some basic preprocessing and feature engineering.","metadata":{}},{"cell_type":"markdown","source":"## A baseline model","metadata":{}},{"cell_type":"code","source":" import os\nfor dirname, _, filenames in os.walk('/kaggle/input'):\n    for filename in filenames:\n        print(os.path.join(dirname, filename))","metadata":{"execution":{"iopub.status.busy":"2021-09-16T12:23:52.939441Z","iopub.execute_input":"2021-09-16T12:23:52.939757Z","iopub.status.idle":"2021-09-16T12:23:52.949115Z","shell.execute_reply.started":"2021-09-16T12:23:52.939715Z","shell.execute_reply":"2021-09-16T12:23:52.948202Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"data=pd.read_csv('/kaggle/input/titanic/train.csv')\ndata['Survived'].value_counts(normalize=True)","metadata":{"execution":{"iopub.status.busy":"2021-09-16T12:23:53.470481Z","iopub.execute_input":"2021-09-16T12:23:53.470764Z","iopub.status.idle":"2021-09-16T12:23:53.501940Z","shell.execute_reply.started":"2021-09-16T12:23:53.470736Z","shell.execute_reply":"2021-09-16T12:23:53.501066Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# preprocess data\ndef cabin_search(row):\n    if row==str(0):\n        return 'no'\n    else:\n        return 'yes'\n\ny = data[\"Survived\"]\nX = data.drop(['Survived','PassengerId'], axis=1)\n\n# add some arbitrary features to make model more complex (XGBoost performs better on complex datasets)\nimport re\n\nX.Cabin.fillna(str(0),inplace=True)\nX['cabin_identified?']=X.apply(lambda x: cabin_search(x['Cabin']),axis=1)\nX['cabin_mark']=data.Cabin.str.extract(r'([A-Z])')\nX['title']=data['Name'].str.extract(r'^[A-Za-z]+,\\s(.+?\\.)')\n\nfrom sklearn.model_selection import train_test_split\n\n# take small train size for demonstration purposes\nX_train_full, X_val_full, y_train, y_val = train_test_split(X, y, train_size=0.6,\n                                                                stratify=y,random_state=0)\n\n# Select categorical columns with relatively low cardinality\ncategorical_cols = [cname for cname in X_train_full.columns if X_train_full[cname].nunique() < 40 and \n                        X_train_full[cname].dtype == \"object\"]\n# Select numerical columns\nnumerical_cols = [cname for cname in X_train_full.columns if X_train_full[cname].dtype in ['int64', 'float64']]","metadata":{"execution":{"iopub.status.busy":"2021-09-16T12:23:54.117488Z","iopub.execute_input":"2021-09-16T12:23:54.117764Z","iopub.status.idle":"2021-09-16T12:23:54.158647Z","shell.execute_reply.started":"2021-09-16T12:23:54.117737Z","shell.execute_reply":"2021-09-16T12:23:54.157720Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Data preprocessing without pipeline (for core XGBoost API) is in the hidden cell below (it's rather long).","metadata":{}},{"cell_type":"code","source":"# preprocess data without pipeline\nfrom sklearn.compose import ColumnTransformer\nfrom sklearn.pipeline import Pipeline\nfrom sklearn.impute import SimpleImputer\nfrom sklearn.preprocessing import OneHotEncoder\n\nnum_imputer=SimpleImputer(strategy='mean')\nX_train_num_imputed=pd.DataFrame(num_imputer.fit_transform(X_train_full[numerical_cols]))\nX_val_num_imputed=pd.DataFrame(num_imputer.transform(X_val_full[numerical_cols]))\n\nX_train_num_imputed.columns=X_train_full[numerical_cols].columns\nX_val_num_imputed.columns=X_val_full[numerical_cols].columns\n\ncat_imputer=SimpleImputer(strategy='most_frequent')\n\nX_train_cat_imputed=pd.DataFrame(cat_imputer.fit_transform(X_train_full[categorical_cols]))\nX_val_cat_imputed=pd.DataFrame(cat_imputer.transform(X_val_full[categorical_cols]))\n\nX_train_cat_imputed.columns=X_train_full[categorical_cols].columns\nX_val_cat_imputed.columns=X_val_full[categorical_cols].columns\n\nOH_encoder=OneHotEncoder(handle_unknown='ignore',sparse=False)\nOH_cols_train = pd.DataFrame(OH_encoder.fit_transform(X_train_cat_imputed))\nOH_cols_val = pd.DataFrame(OH_encoder.transform(X_val_cat_imputed))\n\nOH_cols_train.index = X_train_cat_imputed.index\nOH_cols_val.index = X_val_cat_imputed.index\n\nX_train=pd.concat([X_train_num_imputed,OH_cols_train],axis=1)\nX_val=pd.concat([X_val_num_imputed,OH_cols_val],axis=1)","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2021-09-16T12:23:55.854326Z","iopub.execute_input":"2021-09-16T12:23:55.854754Z","iopub.status.idle":"2021-09-16T12:23:56.157910Z","shell.execute_reply.started":"2021-09-16T12:23:55.854707Z","shell.execute_reply":"2021-09-16T12:23:56.157043Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# specify some hyperparameters for XGBoost, eta stands for learning rate\nparams = {\n            'objective':'binary:logistic',\n            'eta':0.1\n        }   \n\n# transform training and validation set to DMatrix (the main XGBoost dtype)\nX_train.columns=[str(col) for col in X_train.columns]\nX_val.columns=X_train.columns\n\ndmat_train=xgb.DMatrix(X_train,y_train,feature_names=X_train.columns)\ndmat_val=xgb.DMatrix(X_val,y_val,feature_names=X_val.columns)\n\nmodel=xgb.train(params,\n                  dmat_train,\n                  evals=[(dmat_train,'train'),(dmat_val,'validation')])\n\n# transform probabilities to prediction labels (0,1)\ntrain_preds=[1 if pred>0.5 else 0 for pred in model.predict(data=dmat_train)]\nval_preds=[1 if pred>0.5 else 0 for pred in model.predict(data=dmat_val)]\n\nprint('\\nValidation Accuracy : %.2f'%accuracy_score(y_val,val_preds))\nprint('Train Accuracy : %.2f'%accuracy_score(y_train,train_preds))\n\nprint('\\nConfusion Matrix : ')\nprint(confusion_matrix(y_val,val_preds))\n\nprint('\\nClassification Report : ')\nprint(classification_report(y_val,val_preds))","metadata":{"execution":{"iopub.status.busy":"2021-09-16T12:23:57.058851Z","iopub.execute_input":"2021-09-16T12:23:57.059262Z","iopub.status.idle":"2021-09-16T12:23:57.161951Z","shell.execute_reply.started":"2021-09-16T12:23:57.059223Z","shell.execute_reply":"2021-09-16T12:23:57.161255Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plot_confusion_matrix(cm=confusion_matrix(y_val,val_preds),\n                     target_names=['Died','Survived'],\n                     normalize=False)","metadata":{"execution":{"iopub.status.busy":"2021-09-16T12:24:01.693266Z","iopub.execute_input":"2021-09-16T12:24:01.695246Z","iopub.status.idle":"2021-09-16T12:24:01.999349Z","shell.execute_reply.started":"2021-09-16T12:24:01.693706Z","shell.execute_reply":"2021-09-16T12:24:01.998387Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Penalize false negatives with beta = 0.4 (core XGBoost API)\nLet's try to reduce false negatives (37 cases in the baseline model) applying custom loss function with beta = 0.4.\n\nFirst of all, I should define gradients, rules for update custom logloss function and obtain evaluation metric that suits to the function.","metadata":{}},{"cell_type":"code","source":"def first_grad_logreg_beta(predt,dtrain):\n    '''Compute the first derivative for custom logloss function'''\n    y=dtrain.get_label() if isinstance(dtrain,xgb.DMatrix) else dtrain\n    return predt*(beta-beta*y+y)-y\n\ndef second_grad_logreg_beta(predt,dtrain):\n    '''Compute the second derivative for custom logloss function'''\n    y=dtrain.get_label() if isinstance(dtrain,xgb.DMatrix) else dtrain\n    return (y+beta-beta*y)*predt*(1-predt)\n\ndef logregobj_beta(predt,dtrain):\n    '''Custom logloss function update'''\n    grad=first_grad_logreg_beta(predt,dtrain)\n    hess=second_grad_logreg_beta(predt,dtrain)\n    return grad,hess\n\ndef logreg_err_beta(predt,dmat):\n    '''Custom evaluation metric that should be in line with custom loss function'''\n    y=dmat.get_label() if isinstance(dmat,xgb.DMatrix) else dmat\n    predt=np.clip(predt,10e-7,1-10e-7)\n    loss_fn=y*np.log(predt)\n    loss_fp=(1.0-y)*np.log(1.0-predt)\n    return 'logreg_error',np.sum(-(loss_fn+beta*loss_fp))/len(y)","metadata":{"execution":{"iopub.status.busy":"2021-09-16T12:24:17.582329Z","iopub.execute_input":"2021-09-16T12:24:17.583114Z","iopub.status.idle":"2021-09-16T12:24:17.591337Z","shell.execute_reply.started":"2021-09-16T12:24:17.583074Z","shell.execute_reply":"2021-09-16T12:24:17.590173Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"beta=.4\n\n# Dmatrices were taken from the baseline model.\n# new parameters without default objective function and disable default evaluation metric\nparams={'eta':0.1,\n       'disable_default_eval_metric': 1}\n\n# new objective function is obj=logregobj_beta, new metric is feval=logreg_err_beta\nxgb_model=xgb.train(params,\n                   dmat_train,\n                   obj=logregobj_beta,\n                   evals=[(dmat_train,'train'),(dmat_val,'validation')],\n                   feval=logreg_err_beta,\n                    )\n\ntrain_preds=[1 if pred>0.5 else 0 for pred in xgb_model.predict(data=dmat_train)]\nval_preds=[1 if pred>0.5 else 0 for pred in xgb_model.predict(data=dmat_val)]\n\nprint('\\nValidation Accuracy : %.2f'%accuracy_score(y_val,val_preds))\nprint('Train Accuracy : %.2f'%accuracy_score(y_train,train_preds))\n\nprint('\\nConfusion Matrix : ')\nprint(confusion_matrix(y_val,val_preds))\n\nprint('\\nClassification Report : ')\nprint(classification_report(y_val,val_preds))","metadata":{"execution":{"iopub.status.busy":"2021-09-16T12:24:18.199532Z","iopub.execute_input":"2021-09-16T12:24:18.199809Z","iopub.status.idle":"2021-09-16T12:24:18.266628Z","shell.execute_reply.started":"2021-09-16T12:24:18.199782Z","shell.execute_reply":"2021-09-16T12:24:18.266044Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plot_confusion_matrix(cm=confusion_matrix(y_val,val_preds),\n                     target_names=['Died','Survived'],\n                     normalize=False)","metadata":{"execution":{"iopub.status.busy":"2021-09-16T12:24:23.203491Z","iopub.execute_input":"2021-09-16T12:24:23.204000Z","iopub.status.idle":"2021-09-16T12:24:23.489038Z","shell.execute_reply.started":"2021-09-16T12:24:23.203953Z","shell.execute_reply":"2021-09-16T12:24:23.488418Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"So the number of false negative cases reduced from 37 to 27 cases. But overall model performs worse (reduction of  true diseased cases and growth of false positive ones).\n\nLet's take a look on XGBoost sklearn API dealing with false positives.","metadata":{}},{"cell_type":"markdown","source":"## Penalize false positives with beta = 2 (sklearn XGBoost API)\nUser defined functions in sklearn looks different from core XGBoost. Also there is a subtle point relating to a version of sklearn.","metadata":{}},{"cell_type":"code","source":"sklearn.__version__","metadata":{"execution":{"iopub.status.busy":"2021-09-16T12:24:33.434864Z","iopub.execute_input":"2021-09-16T12:24:33.435382Z","iopub.status.idle":"2021-09-16T12:24:33.442358Z","shell.execute_reply.started":"2021-09-16T12:24:33.435348Z","shell.execute_reply":"2021-09-16T12:24:33.441351Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# functions for sklearn 0.23.2. Differences from core XGBoost: order of inputs -  labels go first (except evaluation metric),\n# predt should be transform via logit, because they are inputed as logitraw.\ndef first_grad_logreg_beta_sklearn(dtrain,predt):\n    '''Compute the first derivative for custom logloss function'''\n    y=dtrain.get_label() if isinstance(dtrain,xgb.DMatrix) else dtrain\n    predt=1.0/(1.0+np.exp(-predt))\n    return predt*(beta-beta*y+y)-y\n\ndef second_grad_logreg_beta_sklearn(dtrain,predt):\n    '''Compute the second derivative for custom logloss function'''\n    y=dtrain.get_label() if isinstance(dtrain,xgb.DMatrix) else dtrain\n    predt=1.0/(1.0+np.exp(-predt))\n    return (y+beta-beta*y)*predt*(1-predt)\n\ndef logregobj_beta_sklearn(dtrain,predt):\n    '''Custom logloss function update'''\n    grad=first_grad_logreg_beta_sklearn(dtrain,predt)\n    hess=second_grad_logreg_beta_sklearn(dtrain,predt)\n    return grad,hess\n\ndef logreg_err_beta_sklearn(predt,dmat):\n    '''Custom evaluation metric that should be in line with custom loss function'''\n    y=dmat.get_label() if isinstance(dmat,xgb.DMatrix) else dmat\n    predt=1.0/(1.0+np.exp(-predt))\n    predt=np.clip(predt,10e-7,1-10e-7)\n    loss_fn=y*np.log(predt)\n    loss_fp=(1.0-y)*np.log(1.0-predt)\n    return 'logreg_error',np.sum(-(loss_fn+beta*loss_fp))/len(y)","metadata":{"execution":{"iopub.status.busy":"2021-09-16T12:24:35.172443Z","iopub.execute_input":"2021-09-16T12:24:35.173328Z","iopub.status.idle":"2021-09-16T12:24:35.183381Z","shell.execute_reply.started":"2021-09-16T12:24:35.173281Z","shell.execute_reply":"2021-09-16T12:24:35.182351Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Functions for sklearn 0.22.2.post1 are presented in the hidden cell below (commented).","metadata":{}},{"cell_type":"code","source":"# function for sklearn 0.22.2.post1 (it is google colab version). There is no need to transform predts via logit.\n\n# def first_grad_logreg_sklearn(dtrain,predt):\n#   '''Compute the first derivative for negative logistic regression function'''\n#   y=dtrain.get_label() if isinstance(dtrain,xgb.DMatrix) else dtrain\n#   return predt-y\n\n# def second_grad_logreg_sklearn(dtrain,predt):\n#   '''Compute the second derivative for negative logistic regression function'''\n#   y=dtrain.get_label() if isinstance(dtrain,xgb.DMatrix) else dtrain\n#   return predt*(1-predt)\n\n# def logregobj_sklearn(dtrain,predt):\n#   '''Negative logistic regression function'''\n#   grad=first_grad_logreg_sklearn(dtrain,predt)\n#   hess=second_grad_logreg_sklearn(dtrain,predt)\n#   return grad,hess\n\n# def logreg_err_sklearn(dmat,predt):\n#   y=dmat.get_label() if isinstance(dmat,xgb.DMatrix) else dmat\n#   predt=np.clip(predt,10e-7,1-10e-7)\n#   loss_fn=y*np.log(predt)\n#   loss_fp=(1.0-y)*np.log(1.0-predt)\n#   return 'logreg_error',np.sum(-(loss_fn+beta*loss_fp))/len(y)","metadata":{"_kg_hide-input":true,"execution":{"iopub.status.busy":"2021-09-16T12:24:36.903370Z","iopub.execute_input":"2021-09-16T12:24:36.903994Z","iopub.status.idle":"2021-09-16T12:24:36.908631Z","shell.execute_reply.started":"2021-09-16T12:24:36.903938Z","shell.execute_reply":"2021-09-16T12:24:36.907296Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"beta=2\n\nxgb_classif=xgb.XGBClassifier(objective=logregobj_beta_sklearn,**params)\n\nxgb_classif.fit(X_train,y_train,eval_set=[(X_val,y_val)],\n                eval_metric=logreg_err_beta_sklearn,\n                 early_stopping_rounds=5,\n                verbose=10)\n\ny_pred_val=xgb_classif.predict(X_val)\ny_pred_train=xgb_classif.predict(X_train)\n\nprint('\\nValidation Accuracy : %.2f'%accuracy_score(y_val,y_pred_val))\nprint('Train Accuracy : %.2f'%accuracy_score(y_train,y_pred_train))\n\nprint('\\nConfusion Matrix : ')\nprint(confusion_matrix(y_val,y_pred_val))\n\nprint('\\nClassification Report : ')\nprint(classification_report(y_val,y_pred_val))","metadata":{"execution":{"iopub.status.busy":"2021-09-16T12:24:37.534885Z","iopub.execute_input":"2021-09-16T12:24:37.535515Z","iopub.status.idle":"2021-09-16T12:24:37.646258Z","shell.execute_reply.started":"2021-09-16T12:24:37.535478Z","shell.execute_reply":"2021-09-16T12:24:37.645207Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plot_confusion_matrix(cm=confusion_matrix(y_val,y_pred_val),\n                     target_names=['Died','Survived'],\n                     normalize=False)","metadata":{"execution":{"iopub.status.busy":"2021-09-16T12:24:42.917706Z","iopub.execute_input":"2021-09-16T12:24:42.918040Z","iopub.status.idle":"2021-09-16T12:24:43.223763Z","shell.execute_reply.started":"2021-09-16T12:24:42.918009Z","shell.execute_reply":"2021-09-16T12:24:43.222812Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Count of false postitves cases reduced by nearly half (from 23 to 14 cases). Accuracy remains on the same level.","metadata":{}},{"cell_type":"markdown","source":"# Conclusion\nCustom loss function is a rather powerfull tool that is used to make an impact on FN/FP when standard model is not enough.\nIt may be used with core XGBoost API or Sklearn one.\n\nI prefer the sklearn realisation, because it gives convenient access for pipelines, cross-validation and search for hyperparameters.\n\"Do everything in Sklearn and everyone will be happy\" (c).\n\nWas it an interesting material? Would you like to see LightGBM realisation? What do you prefer core XGBoost or sklearn?\nPlease share your thoughts.","metadata":{}},{"cell_type":"markdown","source":"# References\n1. Overview of binary crossentropy (logloss): https://towardsdatascience.com/understanding-binary-cross-entropy-log-loss-a-visual-explanation-a3ac6025181a\n2. Overview of custom loss function for regression tasks: https://towardsdatascience.com/custom-loss-functions-for-gradient-boosting-f79c1b40466d\n3. XGBoost documentation related to custom loss function: https://xgboost.readthedocs.io/en/latest/tutorials/custom_metric_obj.html\n4. Discussion in stackexchange: https://datascience.stackexchange.com/questions/26972/cost-sensitive-logloss-for-xgboost\n5. Plot_confusion_matrix source: https://www.kaggle.com/grfiv4/plot-a-confusion-matrix\n6. Catch phrase \"Do everything in Sklearn and everyone will be happy\" (c).: https://www.kaggle.com/bextuychiev/write-powerful-code-w-custom-sklearn-transformers","metadata":{}}]}