tutorials/Tutorial_3_Custom_features.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## This tutorial shows how to build custom features in py-boost"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Imports"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os\n",
    "# Optional: set the device to run\n",
    "os.environ[\"CUDA_DEVICE_ORDER\"] = \"PCI_BUS_ID\"\n",
    "os.environ[\"CUDA_VISIBLE_DEVICES\"] = \"0\"\n",
    "\n",
    "os.makedirs('../data', exist_ok=True)\n",
    "\n",
    "import joblib\n",
    "from sklearn.datasets import make_regression\n",
    "\n",
    "# simple case - just one class is used\n",
    "from py_boost import GradientBoosting "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Generate dummy regression data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "CPU times: user 2.34 s, sys: 1.7 s, total: 4.05 s\n",
      "Wall time: 849 ms\n"
     ]
    }
   ],
   "source": [
    "%%time\n",
    "X, y = make_regression(150000, 100, n_targets=10, random_state=42)\n",
    "\n",
    "# we need non negative targets for this example\n",
    "y = y - y.min(axis=0)\n",
    "\n",
    "X_test, y_test = X[:50000], y[:50000]\n",
    "X, y = X[-50000:], y[-50000:]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Custom Loss\n",
    "\n",
    "As it was mentioned in Tutorial_1, not only string alias is valid value for the loss function, but also the instance of Loss class, which is parent class for all loss function\n",
    "\n",
    "#### Now let's build our own MSLE (https://scikit-learn.org/stable/modules/generated/sklearn.metrics.mean_squared_log_error.html) loss function\n",
    "\n",
    "**Note**: Actually we have the built-in MSLE, so you still could use strinng alias for it"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "import cupy as cp\n",
    "from py_boost.gpu.losses import Loss, Metric\n",
    "\n",
    "class CustomRMSLEMetric(Metric):\n",
    "    \"\"\"First, let's define eval metric to estimate model quality while training\"\"\"\n",
    "    \n",
    "    def error(self, y_true, y_pred):\n",
    "        \"\"\"\n",
    "        The simpliest way do define a custom metric is to define .error method\n",
    "        Just tell py_boost how to calculate error at the each point, for out case it is possible\n",
    "        If it is not possible (for ex. ROC-AUC), you should define __call__ method\n",
    "        See the Metric class for the details\n",
    "        \n",
    "        At that stage y_true is already in GPU memory, so we use CuPy to handle it.\n",
    "        Usage is the same as NumPy, just replace np with cp\n",
    "        \n",
    "        Note: the metric is calculated against processed input (see CustomMSLELoss below)\n",
    "        \"\"\"\n",
    "        return (cp.log1p(y_true) - cp.log1p(y_pred)) ** 2\n",
    "    \n",
    "    def compare(self, v0 ,v1):\n",
    "        \"\"\"\n",
    "        The last required method is .compare\n",
    "        It should return True if v0 metric value is better than v1, False othewise\n",
    "        \"\"\"\n",
    "        return v0 < v1\n",
    "    \n",
    "    def __call__(self, y_true, y_pred, sample_weight=None):\n",
    "        \"\"\"\n",
    "        We also update __call__ method to redefine default reduction with square\n",
    "        \"\"\"\n",
    "        return super().__call__(y_true, y_pred, sample_weight) ** .5\n",
    "\n",
    "\n",
    "class CustomMSLELoss(Loss):\n",
    "    \"\"\"Custom MSLE Implementation\"\"\"\n",
    "    \n",
    "    def preprocess_input(self, y_true):\n",
    "        \"\"\"\n",
    "        This method defines, how raw target should be processed before the train starts\n",
    "        We expect y_true has shape (n_samples, n_outputs)\n",
    "        \n",
    "        Here we will not do the actual preprocess, but just check if targets are non negative\n",
    "        \n",
    "        At that stage y_true is already in GPU memory, so we use CuPy to handle it.\n",
    "        Usage is the same as NumPy, just replace np with cp\n",
    "        \n",
    "        Note: All metrics and losses will be computed with this preprocess target\n",
    "        \"\"\"\n",
    "        assert (y_true >= 0).all()\n",
    "        return y_true\n",
    "    \n",
    "    def postprocess_output(self, y_pred):\n",
    "        \"\"\"\n",
    "        Since we modify the target variable, we also need method, that defines \n",
    "        how to process model prediction\n",
    "        \"\"\"\n",
    "        \n",
    "        return cp.expm1(y_pred)\n",
    "    \n",
    "    def get_grad_hess(self, y_true, y_pred):\n",
    "        \"\"\"\n",
    "        This method defines how to calculate gradients and hessians for given loss\n",
    "        Note that training also supports sample_weight, but its applied outside the loss fn,\n",
    "        so we don't need to handle it here\n",
    "        \"\"\" \n",
    "        # grad should have the same shape as y_pred\n",
    "        grad = y_pred - cp.log1p(y_true)\n",
    "        # NOTE: Input could be a matrix in multioutput case!\n",
    "        # But anyway - hessians are ones for all of them\n",
    "        # So, we just create (n_samples, 1) array of ones \n",
    "        # and after that is will be broadcasted over all outputs\n",
    "        # grad should have the same shape as y_pred or (n_samples, 1)\n",
    "        hess = cp.ones((y_true.shape[0], 1), dtype=cp.float32)\n",
    "        \n",
    "        return grad, hess\n",
    "\n",
    "    def base_score(self, y_true):\n",
    "        \"\"\"\n",
    "        One last thing we require to define is base score\n",
    "        This method defines how to initialize an empty ensemble\n",
    "        In simplies case it could be just an array of zeros\n",
    "        But usualy it is better to boost from mean values\n",
    "        Output shape should be (n_outputs, ) \n",
    "        \n",
    "        Note: y_true is already processed array here\n",
    "        \n",
    "        \"\"\"\n",
    "        return cp.log1p(y_true).mean(axis=0)\n",
    "    \n",
    "    \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[14:33:29] Stdout logging level is INFO.\n",
      "[14:33:29] GDBT train starts. Max iter 1000, early stopping rounds 100\n",
      "[14:33:30] Iter 0; Sample 0, score = 0.24603539557907483; \n",
      "[14:33:32] Iter 100; Sample 0, score = 0.1742483282481912; \n",
      "[14:33:35] Iter 200; Sample 0, score = 0.1342659688820449; \n",
      "[14:33:38] Iter 300; Sample 0, score = 0.10731344416487074; \n",
      "[14:33:40] Iter 400; Sample 0, score = 0.08782596307881492; \n",
      "[14:33:43] Iter 500; Sample 0, score = 0.07353079220891415; \n",
      "[14:33:46] Iter 600; Sample 0, score = 0.06300246387723472; \n",
      "[14:33:48] Iter 700; Sample 0, score = 0.05525294291245993; \n",
      "[14:33:51] Iter 800; Sample 0, score = 0.049518312947738174; \n",
      "[14:33:53] Iter 900; Sample 0, score = 0.045306569498698365; \n",
      "[14:33:56] Iter 999; Sample 0, score = 0.04221849120683567; \n",
      "CPU times: user 29.4 s, sys: 4.47 s, total: 33.9 s\n",
      "Wall time: 32.2 s\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<py_boost.gpu.boosting.GradientBoosting at 0x7f9a2070ccd0>"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "%%time\n",
    "model = GradientBoosting(CustomMSLELoss(), CustomRMSLEMetric(), lr=0.01, verbose=100, ntrees=1000)\n",
    "\n",
    "model.fit(X, y, eval_sets=[{'X': X_test, 'y': y_test},])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(50000, 10)"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "model.predict(X_test).shape"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Custom colsample strategy\n",
    "\n",
    "We could also redefine some other things. Let's see the example of creating our bagging strategy. Most of custom things should be done via Callbak. \n",
    "\n",
    "To create callback we should inherit Callbak class. There are 4 methods, that could be redefined:\n",
    "        - before_train - outputs None\n",
    "        - before_iteration - outputs None\n",
    "        - after_train - outputs None\n",
    "        - after_iteration - outputs bool - if training should be stopped after iteration\n",
    "\n",
    "    Methods receive build_info - the state dict, that could be accessed and modifier\n",
    "\n",
    "    Basic build info structure:\n",
    "\n",
    "    build_info = {\n",
    "            'data': {\n",
    "                'train': {\n",
    "                    'features_cpu': np.ndarray - raw feature matrix,\n",
    "                    'features_gpu': cp.ndarray - uint8 quantized feature matrix on GPU,\n",
    "                    'target': y - cp.ndarray - processed target variable on GPU,\n",
    "                    'sample_weight': cp.ndarray - processed sample_weight on GPU or None,\n",
    "                    'ensemble': cp.ndarray - current model prediction (with no postprocessing,\n",
    "                        ex. before sigmoid for logloss) on GPU,\n",
    "                    'grad': cp.ndarray of gradients on GPU, before first iteration - None,\n",
    "                    'hess': cp.ndarray of hessians on GPU, before first iteration - None,\n",
    "\n",
    "                    'last_tree': {\n",
    "                        'leaves': cp.ndarray - nodes indices of the last trained tree,\n",
    "                        'preds': cp.ndarray - predictions of the last trained tree,\n",
    "                    }\n",
    "\n",
    "                },\n",
    "                'valid': {\n",
    "                    'features_cpu' the same as train, but list, each element corresponds each validation sample,\n",
    "                    'features_gpu': ...,\n",
    "                    'target': ...,\n",
    "                    'sample_weight': ...,\n",
    "                    'ensemble': ...,\n",
    "\n",
    "                    'last_tree': {\n",
    "                        'leaves': ...,\n",
    "                        'preds': ...,\n",
    "                    }\n",
    "\n",
    "                }\n",
    "            },\n",
    "            'borders': list of np.ndarray - list or quantization borders,\n",
    "            'model': GradientBoosting - model, that is trained,\n",
    "            'mempool': cp.cuda.MemoryPool - memory pool used for train, could be used to clean memory to prevent OOM,\n",
    "            'builder': DepthwiseTreeBuilder - the instance of tree builder, contains training params,\n",
    "\n",
    "            'num_iter': int, current number of iteration,\n",
    "            'iter_scores': list of float - list of metric values for all validation sets for the last iteration,\n",
    "        }\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "import cupy as cp\n",
    "from py_boost.callbacks.callback import Callback\n",
    "\n",
    "class ColumnImportanceSampler(Callback):\n",
    "    \"\"\"\n",
    "    This class implements a sampling strategy, \n",
    "    that sample columns in proportion to thier importance at each step\n",
    "    \n",
    "    We should implement __call__ method to use it as sampler\n",
    "    \"\"\"\n",
    "    def __init__(self, rate=0.5, smooth=0.1, \n",
    "                 update_freq=10, inverse=False):\n",
    "        \"\"\"\n",
    "        \n",
    "        Args:\n",
    "            rate: float, sampling rate\n",
    "            smooth: float, smoothing parameter\n",
    "            update_freq: int importance update frequency\n",
    "            inverse: inverse the probability of sampling\n",
    "\n",
    "        Returns:\n",
    "\n",
    "        \"\"\"\n",
    "        # Custom columnns sampler based on feature importance\n",
    "        self.rate = rate\n",
    "        self.smooth = smooth\n",
    "        self.update_freq = update_freq\n",
    "        self.inverse = inverse\n",
    "        \n",
    "    def before_iteration(self, build_info):\n",
    "        \"\"\"\n",
    "        Define what should be doe before each iteration\n",
    "        \"\"\"\n",
    "        # Update feature importance\n",
    "        num_iter = build_info['num_iter']\n",
    "        \n",
    "        if (num_iter % self.update_freq) == 0:\n",
    "            # update probabilities with actual importance\n",
    "            p = build_info['model'].get_feature_importance() + 1e-3\n",
    "            p = cp.asarray(p) / (p.sum())\n",
    "            # inverse if needed\n",
    "            if self.inverse:\n",
    "                p = 1 - p\n",
    "                p = p / p.sum()\n",
    "            # apply smoothing\n",
    "            self.p = p * (1 - self.smooth) + cp.ones_like(p) * self.smooth / p.shape[0]\n",
    "            \n",
    "    def __call__(self):\n",
    "        \"\"\"\n",
    "        Method should return the array of indices, that will be used\n",
    "        to grow the tree at the current step\n",
    "        \"\"\"\n",
    "        # Sample rows\n",
    "        n = self.p.shape[0]\n",
    "        index = cp.random.choice(cp.arange(n, dtype=cp.uint64), \n",
    "            size=int(self.rate * n), p=self.p)\n",
    "        \n",
    "        return index"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[14:33:57] Stdout logging level is INFO.\n",
      "[14:33:57] GDBT train starts. Max iter 1000, early stopping rounds 100\n",
      "[14:33:57] Iter 0; Sample 0, score = 0.24644921389665553; \n",
      "[14:33:59] Iter 100; Sample 0, score = 0.17590711477798346; \n",
      "[14:34:00] Iter 200; Sample 0, score = 0.13484779001390923; \n",
      "[14:34:02] Iter 300; Sample 0, score = 0.10826939489014992; \n",
      "[14:34:03] Iter 400; Sample 0, score = 0.08943576705705947; \n",
      "[14:34:05] Iter 500; Sample 0, score = 0.0753772653073726; \n",
      "[14:34:07] Iter 600; Sample 0, score = 0.06446689810058637; \n",
      "[14:34:08] Iter 700; Sample 0, score = 0.05593631183289121; \n",
      "[14:34:10] Iter 800; Sample 0, score = 0.04973719737069171; \n",
      "[14:34:12] Iter 900; Sample 0, score = 0.045194617065396514; \n",
      "[14:34:13] Iter 999; Sample 0, score = 0.041950662857688406; \n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<py_boost.gpu.boosting.GradientBoosting at 0x7f9a2070ca30>"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# create model with new sampler   \n",
    "# if we pass new sampler to the colsample argument it will used instead of default\n",
    "# it will also be added to the callback pipeline automatically\n",
    "# you should not pass samplers to the callbacks argument\n",
    "\n",
    "model = GradientBoosting(CustomMSLELoss(), CustomRMSLEMetric(), \n",
    "                         colsample=ColumnImportanceSampler(0.5), \n",
    "                         lr=0.01, verbose=100, ntrees=1000 )\n",
    "\n",
    "model.fit(X, y, eval_sets=[{'X': X_test, 'y': y_test},])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(50000, 10)"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "model.predict(X_test).shape"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "rapids-22.06",
   "language": "python",
   "name": "rapids-22.06"
  },
  "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.9.13"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}