Mismatch is not uncommon
This 'mismatch' occurs, for example, in any case where GLM is being tuned via cross validation (an example is the R function cv.glm
).
- The GLM model is fitted by maximizing the likelihood function.
- The nuisance parameters of the model may be tuned by a different function, e.g. $R^2$, or some other measure of performance.
The mismatch occurs whenever the objective for the fitting of coefficients (e.g. likelihood) is different from the objective for model selection (e.g. predictive power or some specialized loss function).
Why two different loss functions?
The objective for fitting of the coefficients follows the statistical model for the generation of the data and is done in order to filter the model from the noise (the first cost function weighs the noise). The objective for selecting the model may be different and is done in order to optimize some performance measure (the second cost function weights the bias).
Is the fitted model gonna perform better when we make a fit that is closer to the true data generation model (first loss function), or a fit that is closer to some performance measure (second loss function)?
Whether it is good or bad to have a mismatch between the two loss functions will depend on the relative importance of the deterministic part and the random part of the model.
- The deterministic part: How close is the model to the true relationship (the amount of bias)?
- The random part: How much noise is present?
Using a different loss function for the fitting than for the model selection will emphasize a better filtering of the model from the noise, this may be at the cost of less reduction in the bias in the model (as far as this bias can be reduced by choosing different coefficients).
Reversal of the question
The better question is:
Would it be justified to use the model selection (testing) objective in fitting the model coefficients (training)?
Do we get better performance when we fit the model in training to minimize the performance objective instead of some more 'natural' way to fit the model?
For instance when we have data generated according to some Poisson process then it would be 'natural' to fit it using Poisson regression. But, should we fit it with a least squares fit instead when we are (for whatever reason) measuring performance by the least squares distance?
If no or little bias is expected (ie little discrepancy between true relationship of the mean and modeled for the relationship of the mean) then it would be better to fit a Poisson regression, if you like to minimize the least squares of a prediction.
###Two computational examples
The example below uses glmnet
in R to compare a Poisson regression with a least squares regression.
The least squares regression 'loses' 65 times out of 100 in getting the lowest least squares error.
library(glmnet)
set.seed(1)
n=1000
p=100
nzc=trunc(p/10)
# keeping count how often the one method performs better
Poisson_vs_Gauss = 0
for (i in 1:100) {
# make random matrix X and independent variable Y
x=matrix(rnorm(n*p),n,p)
beta=rnorm(nzc)
fx= x[,seq(nzc)] %*% beta
mu=exp(fx/10)
y=rpois(n,mu)
# perform penalized GLM in two ways
cvob1=cv.glmnet(x,y,type.measure="mse")
cvob4=cv.glmnet(x,y,type.measure="mse",family="poisson")
# compare the two
if (min(cvob1$cvm)>min(cvob4$cvm)) {
Poisson_vs_Gauss = Poisson_vs_Gauss+1
}
}
plot(log(cvob1$lambda),cvob1$cvm,
ylim=c(1,1.2),
xlab = "log(lambda)", ylab = "MSE",
pch=21,col="black",bg="white",cex=0.7)
points(log(cvob4$lambda),cvob4$cvm,
pch=21,col="black",bg="gray",cex=0.7)
A problem with this example is the link function being different. I could not find standard libraries that allow computation of penalized glm with adjusted link functions.
Therefore the following example where the polynomial degree is adjusted to model some function.
The least squares regression 'loses' 64 times out of 100 in getting the lowest least squares error.
library(boot)
set.seed(1)
#data
# keeping count how often the one method performs better
Poisson_vs_Gauss_2 = 0
for (i in 1:100) {
# create data
# X is a polynomial of order 1 to 5
# Y is Poisson distributed data with a mean modeled by a sinus function.
x <- seq(1,10,length.out=40)
y <- rpois(length(x),5+4*sin(x/10*pi))
datat <- data.frame(y=y,
x0 = rep(1,length(x)),
x1 = x,
x2 = x^2,
x3 = x^3,
x4 = x^4,
x5 = x^5)
# computing performance for two models
rms1 = rep(NA, 5)
rms2 = rep(NA, 5) # init result vectors
degree = 1:5
for (d in degree) {
# Gaussian model (minimizes squared error)
fit1 <- glm(y ~ 0 + .,
data = datat[,1:(2+d)],
family = gaussian(link="identity"))
rms1[d] <- cv.glm(data = datat[,1:(2+d)], glmfit = fit1)$delta[1]
# Poisson model (minimizes likelihood)
fit2 <- glm(y ~ 0 + .,
data = datat[,1:(2+d)], start = fit1$coefficients,
family = poisson(link="identity"))
rms2[d] <- cv.glm(data = datat[,1:(2+d)], glmfit = fit2)$delta[1]
}
if (min(rms1)>min(rms2)) {
Poisson_vs_Gauss_2 = Poisson_vs_Gauss_2+1
}
}
Poisson_vs_Gauss_2
plot(degree, rms1,
ylim=c(min(rms1,rms2),max(rms1,rms2)),
xlab = "log(lambda)", ylab = "MSE",
pch=21,col="black",bg="white",cex=0.7)
points(degree, rms2,
pch=21,col="black",bg="gray",cex=0.7)