Quantile Loss Function Implementation

Standard regression techniques often focus on predicting the conditional mean of the target variable. Quantile regression, however, allows for modeling conditional quantiles (like the median, 10th percentile, or 90th percentile), providing a more comprehensive understanding of the relationship between predictors and the target variable's distribution. To perform quantile regression using gradient boosting, it is necessary to define and implement the appropriate loss function: the quantile loss, also known as the pinball loss.

The Quantile Loss Function

The quantile loss function for a specific quantile $\alpha \in (0, 1)$ measures the error between the true value $y$ and the predicted quantile $\hat{y}$. It is defined as:

$$L_\alpha(y, \hat{y}) = \begin{cases} \alpha (y - \hat{y}) & \text{if } y - \hat{y} > 0 \\ (1 - \alpha) (\hat{y} - y) & \text{if } y - \hat{y} \le 0 \end{cases}$$

This can be written more compactly using the indicator function $I(\cdot)$:

$$L_\alpha(y, \hat{y}) = (y - \hat{y}) (\alpha - I(y - \hat{y} < 0))$$

Notice the asymmetry:

• When $\alpha = 0.5$ (the median), the loss simplifies to $0.5 |y - \hat{y}|$, which is proportional to the Mean Absolute Error (MAE). The penalty for underprediction ($\hat{y} < y$) is equal to the penalty for overprediction ($\hat{y} > y$).
• When $\alpha > 0.5$, the loss penalizes underpredictions ($y > \hat{y}$) more heavily than overpredictions. This encourages the model to predict higher values, corresponding to upper quantiles.
• When $\alpha < 0.5$, the loss penalizes overpredictions ($y < \hat{y}$) more heavily, pushing predictions towards lower quantiles.

The following plot visualizes the pinball loss for different values of $\alpha$.

The plot shows the characteristic 'pinball' shape of the quantile loss. Note how the slope changes at zero error, and the asymmetry depends on the quantile $\alpha$. For $\alpha=0.9$, positive errors (underpredictions) are penalized more, while for $\alpha=0.1$, negative errors (overpredictions) incur a larger penalty. $\alpha=0.5$ treats positive and negative errors symmetrically.

Gradients and Hessians for Boosting

Gradient boosting algorithms iteratively fit base learners (typically trees) to the negative gradient of the loss function with respect to the current prediction. Therefore, we need the first derivative (gradient) of $L_\alpha(y, \hat{y})$ with respect to $\hat{y}$. Some advanced implementations like XGBoost and LightGBM also use the second derivative (Hessian) for faster convergence and regularization.

Let's find the gradient:

$$g = \frac{\partial L_\alpha(y, \hat{y})}{\partial \hat{y}} = I(y - \hat{y} < 0) - \alpha = I(\hat{y} > y) - \alpha$$

So, the gradient is:

• $-\alpha$ if $\hat{y} \le y$ (underprediction or correct prediction)
• $1 - \alpha$ if $\hat{y} > y$ (overprediction)

Now, consider the Hessian:

$$h = \frac{\partial^2 L_\alpha(y, \hat{y})}{\partial \hat{y}^2} = \frac{\partial g}{\partial \hat{y}}$$

The gradient is a step function, which means its derivative is zero everywhere except at the point $\hat{y} = y$, where it is undefined (mathematically, it involves a Dirac delta function). Standard gradient boosting implementations, particularly those using second-order approximations like XGBoost, require a defined Hessian.

How do we handle this?

1. Built-in Objectives: Libraries like LightGBM often have built-in support for quantile regression (objective='quantile'). In these cases, the library handles the Hessian internally, possibly using approximations or specific algorithms suited for non-smooth objectives. Using the built-in objective is generally the recommended approach if available.
2. Custom Objective - Approximate Hessian: When defining a custom objective function, you must return both gradient and Hessian values for each data point. Since the true Hessian is problematic, a common practice is to return a small positive constant (e.g., 1, or even smaller values like 1e-6) for the Hessian. This provides numerical stability and allows the algorithm (which often uses $g/h$ or similar terms) to proceed. While not strictly the true second derivative, it often works adequately in practice.
3. Custom Objective - Zero Hessian: Some frameworks might allow or implicitly handle a zero Hessian, effectively falling back to first-order gradient methods for those steps. However, explicitly providing a small positive value is often more effective in libraries expecting second-order information.

Implementing Custom Quantile Objectives

Let's see how you would structure a custom quantile loss function for libraries like XGBoost or LightGBM in Python. These libraries expect a function that takes the current predictions (preds) and the true labels (dtrain, which contains y) and returns the gradient and Hessian for each sample.

Here's a Python function structure for a custom quantile objective:

import numpy as np

def quantile_objective(alpha):
    """
    Custom objective function for quantile regression.

    Parameters:
    alpha (float): The target quantile, must be in (0, 1).

    Returns:
    callable: A function compatible with XGBoost/LightGBM custom objectives.
    """
    def objective_function(preds, dtrain):
        """
        Calculates gradient and Hessian for quantile loss.

        Parameters:
        preds (np.ndarray): Current model predictions.
        dtrain: Data container (e.g., xgboost.DMatrix or lightgbm.Dataset)
                containing true labels via dtrain.get_label().

        Returns:
        grad (np.ndarray): The gradient of the loss with respect to preds.
        hess (np.ndarray): The Hessian of the loss with respect to preds.
        """
        labels = dtrain.get_label()
        errors = preds - labels # Note: Using preds - labels aligns with I(preds > y) convention

        # Calculate Gradient
        grad = np.where(errors > 0, 1 - alpha, -alpha)

        # Calculate Hessian (using a small constant approximation)
        # A small positive constant helps ensure stability in the algorithm.
        # The exact value might require tuning or experimentation.
        hess = np.full_like(preds, 1.0) # Or a smaller value like 1e-6

        return grad, hess
    return objective_function

# Example Usage
# alpha_value = 0.75 # Target the 75th percentile
# custom_obj = quantile_objective(alpha_value)

# In XGBoost:
# model = xgb.train(params, dtrain, num_boost_round=100, obj=custom_obj)

# In LightGBM:
# model = lgb.train(params, dtrain, num_boost_round=100, fobj=custom_obj)

Important Considerations:

• Library Support: Always check the documentation of your chosen library (XGBoost, LightGBM, CatBoost). They might offer built-in quantile regression objectives (e.g., objective='quantile' with alpha parameter in LightGBM) which are often optimized and preferred over custom implementations. CatBoost also supports quantile loss functions.
• Hessian Approximation: The choice of the constant value for the Hessian in a custom objective can sometimes influence convergence speed and stability. If using a custom implementation, you might experiment with different small positive values (e.g., 1.0, 0.1, 1e-3, 1e-6).
• Tuning $\alpha$: The quantile level $\alpha$ is treated like a hyperparameter if you are trying to model a specific quantile based on domain knowledge or problem requirements. If you need to model multiple quantiles, you typically train separate models for each desired $\alpha$.
• Evaluation Metrics: When training a quantile regression model, standard metrics like RMSE (which focuses on the mean) are not appropriate for evaluation. Use metrics relevant to quantiles, such as the quantile loss itself on a validation set, or metrics like the Winkler score if predicting intervals. You might need to implement custom evaluation metrics as well (covered in Chapter 7).

By implementing or utilizing the quantile loss function, you equip gradient boosting models with the capability to move past mean prediction and estimate conditional quantiles, offering much richer insights into the underlying data distribution for various specialized tasks.