Linear Regression in Spark's MLLib gives a seeming incorrect result

I am running the example found here.

The training data for the model can be found in this CSV, where the first column is the response variable and the second column is a space separated list of predictors.

After running the example with the modified parameters of numIterations = 1000 and stepSize = 1, the result it gives is:

intercept = 0.0
weights   = [0.5047290722337779, 0.1920796185138009, 0.3055517342539938, 0.10461255324752054, 0.4918056810542394, -0.21803087595354512, -0.059099512695448504, 0.0]


The trained model is very poor. Furthermore, it is different from the exact solution given by (XTX)-1XTy for least-squares problem. The actual solution using this formula is:

intercept = 2.46493292212375
weights   = [0.679528141237974, 0.263053065732544, -0.141464833536172, 0.210146557221827, 0.305200597125096, -0.288492772453545, -0.021305038802947, 0.266955762119923]


Nothing I've seen in the documentation of LinearRegressionWithSGD says the model is not an OLS model, so I assume it is. The documentation here says the model minimizes the MSE.

The model has a null intercept and the last feature is dropped from the model, which makes me think Spark is internally normalizing the data, training the model on the normalized data, and returning that model. Does anyone have any experience with this?

I have tried comparing this result those provied by lm in R and the linear regression in the Spark's ML library. Both of them return the latter result.

> lm(y ~ x, data)

Call:
lm(formula = y ~ x, data = data)

Coefficients:
(Intercept)          xV2          xV3          xV4          xV5          xV6          xV7          xV8          xV9
2.46493      0.67953      0.26305     -0.14146      0.21015      0.30520     -0.28849     -0.02131      0.26696


If I force the fitting in R to also lack an intercept term, the result is still different:

Call:
lm(formula = y ~ x - 1, data = data)

Coefficients:
xV2      xV3      xV4      xV5      xV6      xV7      xV8      xV9
0.5999   0.1859   0.2808   0.1108   0.4003  -0.5932  -0.6133   0.9169

• The solution provided by $(X^T X)^{-1} X^T y$ is not always numerically stable. I believe most statistical software packages use a Cholesky decomposition. I recommend running a regression in R or SPSS and comparing those results with what Spark provides. – Jon Feb 1 '17 at 18:33
• Running lm on the raw data in R produces the latter result. I'll edit my post to include more details. – Jon Claus Feb 1 '17 at 19:23
• Could you include the results from R as comparison? – Jon Feb 1 '17 at 19:27
• Can you please update the R output to include a model described by hxd1011? lm(formula = y ~ x-1, data = data). Otherwise, your code compares two different models. – Jon Feb 1 '17 at 19:57
• @Jon FYI R's lm uses QR decomposition by default. – Sycorax Feb 1 '17 at 21:40

1. In your spark results, the intercept is $0$, are you enforcing the model pass though the origin? (for getting to get all $1$ column in your data)? try this in R lm(formula = y ~ x-1, data = data) to see what coefficient you get.