Ensemble model of simple univariate linear regression models

Question

Calculating univariate linear regression and correlation is simple in SQL.

It can be done like so:

select 
avg(y) - avg(x) *
((count(*) * sum(x*y)) - (sum(x)*sum(y)))/
((count(*) * sum(x*x)) - (Sum(x)*Sum(x)))   as b0, --intercept
((count(*) * sum(x*y)) - (sum(x)*sum(y)))/
((count(*) * sum(x*x)) - (sum(x)*sum(x))) AS b1,
(Avg(x * y) - Avg(x) * Avg(y)) / (StDevP(x) * StDevP(y)) AS Correlation
from SomeTable

In my scenario, we already have a nominal variable with about 100 values. However, in each of those groups we would like to differentiate the population even further.

Let's say one has a table of 500 numeric continuous potential predictor variables and one target response variable.

For business reasons, we would like to keep using that nominal category variable but apply a simple linear regression after that. Of the 500 additional variables, there are ones that are good and most that are junk for predicting the target variable. From doing regression analysis on a sample, I know which ones overall that lasso would pick in addition to the nominal category. For simplicity, we would like to process this in SQL rather than have a job put chunks through R, Python etc.

What is the best way to choose which simple linear regression models to compose in the ensemble model? (The method of combining being taking the average of all chosen univariate regression models.)

Correlation feature selection chooses predictor variables that are highly correlated with the classification, yet uncorrelated to each other.

So one could pick the linear models based to be included in the ensemble by choosing those models in order of best correlation which are not correlated with previously chosen models. However, this has the complexity of computing a correlation matrix.

More simply one could do a step-wise forward selection and chose those models which improve the ensemble average.

Right now this is just a thought experiment, but I am wondering if someone else has researched this more fully?

Are there any systems that work like this?

What are the massive short-comings?

Obviously, Spark and the like are better, but under considered is a plain SQL solution.

References:

https://stackoverflow.com/questions/2536895/are-there-any-linear-regression-function-in-sql-server/3283119#3283119

http://www.experts-exchange.com/questions/26516795/Calculating-correlation-coefficient-in-SQL-Server.html

Answer 1

I would suggest an alternative if you insist on doing simple univariate regression... Regress on prediction errors. ie predict using one variable, now do univariate regression of error on next variable... This is equivalent to multiple regression but not so numerically well behaved.

There is another alternative Elements of statistical learning, 3.2.3 Multiple regression from simple univariate regression (p.52), orthogonalise the independent variables rather than the dependent variable. I am afraid I haven't got my head round whether they are equivalent or not!

Answer 2

I tested @seanv507 suggestion and it works in R code as a poc. We didn't want to shuttle a large data sets around, so this will work in SQL.

options(scipen=10)
inc <- function(x,y=1){
    if(is.null(x)) {
        eval.parent(substitute(x <- y))
    }else{
        eval.parent(substitute(x <- x + y))
    }
}
rmse=function(x,y,k=0){
  return( sqrt(sum((x-y)^2)/(length(x)-k)));
}
yvar='y'
df=data.frame(x1=runif(1000*10),x2=rlnorm(1000),x3=rnorm(1000),x4=rnorm(1000))
df$y=with(df,x1+2*x2+4*x3)+runif(1000)

cdf=as.data.frame(cor(df,df[[yvar]]))
result=list('(Intercept)'=0)
df$Temp=df[[yvar]]
    for(n in row.names(cdf[order(-abs(cdf$V1)),,drop=F])){
    if(n==yvar) next
    m=lm(as.formula(paste('Temp ~',n)),df)
    inc(result[['(Intercept)']],coef(m)[1])
    inc(result[[n]],coef(m)[n])
    df$Temp=m$residuals
}
df$Temp=NULL
    m=lm(y~.,df)
    summary(m)
    for(n in names(coef(m))){
        message(sprintf("%-15s %f",n,result[[n]]))
    }
    df$fit=predict(m)
for(n in names(coef(m))){
    if(n=='(Intercept)') {
        df$fit2=rep(result[[n]],nrow(df))
        }else{
            df$fit2=df$fit2+df[[n]]*result[[n]]
        }
    }
    plot(df$fit,df$fit2)
    grid()
    rmse(df$fit,df$fit2)
    rmse(df$fit,df$fit2)/mean(df$fit) #RSD AKA CV

Running, I get nearly the same coef. It is not perfect, but it should work for my scenario.

Call:
lm(formula = y ~ ., data = df)

Residuals:
     Min       1Q   Median       3Q      Max 
-0.50973 -0.25359  0.00855  0.25565  0.51839 

Coefficients:
             Estimate Std. Error  t value Pr(>|t|)    
(Intercept) 0.4914712  0.0062200   79.015   <2e-16 ***
x1          1.0103055  0.0102303   98.756   <2e-16 ***
x2          2.0003738  0.0011939 1675.451   <2e-16 ***
x3          4.0089037  0.0029015 1381.644   <2e-16 ***
x4          0.0008583  0.0028330    0.303    0.762    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.2938 on 9995 degrees of freedom
Multiple R-squared:  0.9978,    Adjusted R-squared:  0.9978 
F-statistic: 1.137e+06 on 4 and 9995 DF,  p-value: < 2.2e-16
(Intercept)     0.602780
x1              1.018233
x2              1.931142
x3              4.007331
x4              -0.002166
Out[8]:
0.170448782159297
Out[8]:
0.0403919453688196