Trjtdtk

I am looking at a software development pipeline where I am predicting the lead time of different products flowing through the pipeline.

After applying a boxcox transformation on the lead time (target variable) and creating a XGBoost regressor model I can see that the residuals follow a t-locationScale distribution.

So now I looked at this guide which describes a method to create a prediction interval for any regression model assuming that the residuals are normally distributed. https://qucit.com/a-simple-technique-to-estimate-prediction-intervals-for-any-regression-model_en/

But I tried to tweak it to my distribution.

So a t-locationScale distribution has a $sigma$, $mu$ and $nu$ parameter. The variance is only defined for $nu>2$. My specific distribution has $nu = 2.56$ and $mu = 0.04$, $sigma = 0.97$ So I could take the 95% interval of this distribution and say that for any $haty$, the prediction interval is the 95% interval of the residual distribution.

But I want to take into consideration that the prediction interval should change with different inputs. I created a regressor model, which I trained and then made predictions using the validation set. I then took the square of the error and trained an additional error model on this data. Such that the error model could predict the variance of the residuals distribution.

 xgb = XGBoostRegressor()
 xgb.fit(X_train,y_train)
 y_hat = xgb.predict(X_val)
 val_error = (y_hat-y_val)**2

 xgb_error = XGBoostRegressor()
 xgb_error.fit(X_val, val_error)

 variance_hat_residuals = xgb_error.predict(X_test)

The relationship between variance and $sigma$ and $nu$ for a t-locationScale distribution is

var = $sigma^2 *fracnunu-2$

Now here is where I make an assumption which I am not sure makes sense.

I assume that the degrees of freedom $nu$ is the same as for all residuals, $nu = 2.56$ and then I solve for $sigma$ through the following.

$hatsigma = sqrtfrachatvar*(nu-2)nu$

And estimate the lower and upper quantiles from this distribution.

 residual_distribution_lower_quantile = scipy.stats.t.ppf(q = 0.025, df = 2.56, scale = sigma)
 residual_distribution_upper_quantile = scipy.stats.t.ppf(q = 0.0975, df = 2.56, scale = sigma)

I then predict the lead time $haty$ and say that the mean of the distribution is $haty$

 pred = xgb.prediction(X_test)
 lower_interval = pred + residual_distribution_lower_quantile
 upper_interval = pred + residual_distribution_upper_quantile

Does it make sense to make the claim of $nu$ is static? My score for the prediction interval is now $81%$ since I am clearly simplifying the problem.

Any suggestions for improving my method?