# Bias correction for regression with t-distributed error

I have a GAM /regression model which is originally defined as:

log10(Y)~s(log10(X1))+s(log10(X2))+s(log10(X3)) #using R mgcv


The response needs to be back transformed. As simply doing exp(Y) would give a biased result, I do the bias correction via: Y_unbiased~exp(Y)*exp(((log(10)*sd(model))^2)/2) see: https://derekogle.com/FSA/reference/logbtcf.html

Now it seems that for the data on log10-log10 the assumption of a gaussian error is not correct and I still have very heavy tails. The scat() family from mgcv (https://rdrr.io/cran/mgcv/man/scat.html) provides optimal qqplots for the model. It is based on the t-distribution.

How would now a bias correction look like for a t-distribution ? Is there one at all?

P.S. I also tried to model:

Y~s(log10(X1))+s(log10(X2))+s(log10(X3)) #using R gamlss and a lognormal error distribution
Y~s(log10(X1))+s(log10(X2))+s(log10(X3)) #using R mgcv and a Gamma error distribution with log link


none of those gave as good models i.e. qqplots as the scat family on log10(Y) so how would a bias correction for scat() i.e. the t-distribution look like?

• What kind of data do you have? Is it a time series? Or is it panel data? Mar 25 at 18:14
• it is a time series, but without any autocorrelation. but this should be more general question, how to backtransform for the case scat() heavy taled data works better than gaussian on a log(Y) fit. Mar 25 at 19:50
• The expectation $E[e^X]$ can be found using the moment generating function. For the t-distribution this function does not exist. Also computing directly with Wolfram Alpha gives no result. It is in this term $(1+x/y)^{-y}$ which approaches $exp(-x)$ for $y \to \infty$, but if you go far enough into the tail, then it does not approach zero as fast as exponentiation and $exp(x)$ will not overtake $exp(-x^2)$ but it does overtake $(1+x^2/y)^{-y}$. So $E[e^X]=\infty$. Mar 26 at 7:53
• Thus a mean-unbiased estimator can not be defined because the distribution of $Y=exp(X)$ has infinite mean when $X$ is t-distributed. You can still speak about median-biasedness though. Mar 26 at 8:09
• thanks a lot! About the median-biasedness: Is exp(Y) still the unbiased median estimate of the distribution? and/or can be seen as the geometric mean? comparable to exp(mu) in the lognormal distribution Mar 26 at 23:26

The inherent issue with your model is that your variables are not fully differentiated. Therefore, they all have unit roots, are not stationary, and lead to overall mispecified models with residuals that probably breach all the assumptions of typical models including residuals not being stationary, being heteroskedastic, autocorrelated, not normally distributed, etc.

The solution is to transform your variables. First, instead of Log, take the first difference in such logs from one period to another. This should solve most of the mentioned problems.

If you want to further stabilize the data, you could make additional variable transformations such as:

1. Standardize the data, so each variable is on a comparable scale or # of standard deviations away from the average of that same variable;

2. Normalize the data, so each observed value ranges between 0 (minimum value) and 1 (maximum value) or use another normalization where - 1 (minimum value);

3. Remove the seasonality effect from the data.

4. Remove the heteroskedasticity from the data if the standard deviation (or technically the variance) changes a lot over time.

For an excellent short video that shows how to do a lot of those transformations in Python, check the URL below https://www.youtube.com/watch?v=7_Js8h709Dw

Here is an equivalent video for variable transformation using R. https://www.youtube.com/watch?v=sWDYqnh376o

It is very likely that once you improve the structure of your variables and explore the optimal variable transformation, all your problems related to model biases will go away.

• thanks for the answer. I agree that rescaling may help , but when e.g. z-scoring I always have to remember mean and sd to transform new measurements of my covariates. what if the sd or mean of my covariates change over time e.g. sd(olddata) is significantly different from sd(newdata) or sd(olddata,newdata) , does this make thinks even more complicated? Mar 26 at 23:56
• update: I got nice results when using gaulss() as error family fitting both mean and sd dependent on my covariates. Inserting that into exp(Y)* exp(((log(10)*sd(model))^2)/2) where sd is now not fix but dependent on my covariates gives reasonable results for my data. Mar 27 at 1:15