Neural net without hidden layers should be a simple linear model: why do I get so different results? I’m reading „Elements of Statistical Learning“ where Hastie et al. describe in Section 11.3 on neural nets (p. 394), that (in short) if there are no hidden layers in a neural net (so without non-linear transformation of features), then the model collapses to a linear model.

I tested this in a Keras model with only the output layer included in the model. I noticed that the results (mean absolute error, mae) differ a lot contingent on the choice of the learning rate and the number of training epochs. For some combinations of training epochs and learning rate, I get a mae which is very close to the one of an ordinary least square (OLS) model. 

However, in many cases, results are hugely different (and of course worse than under OLS). 
I don’t understand why this happens. Can someone point me to an answer?

Here is my playcode (using Boston housing data):

library(keras)
library(tibble)
# Data
boston_housing &lt;- dataset_boston_housing()

c(train_data, train_labels) %&lt;-% boston_housing$train
c(test_data, test_labels) %&lt;-% boston_housing$test
paste0("Training entries: ", length(train_data), ", labels: ", length(train_labels))
train_data[1, ] # Display sample features, notice the different scales

column_names &lt;- c('CRIM', 'ZN', 'INDUS', 'CHAS', 'NOX', 'RM', 'AGE', 
                  'DIS', 'RAD', 'TAX', 'PTRATIO', 'B', 'LSTAT')
train_df &lt;- as_tibble(train_data)
colnames(train_df) &lt;- column_names
train_df
train_labels[1:10] # Display first 10 entries


############################################
# Feature engineering

# Test data is *not* used when calculating the mean and std.
# Normalize training data
train_data &lt;- scale(train_data) 

# Use means and standard deviations from training set to normalize test set
col_means_train &lt;- attr(train_data, "scaled:center") 
col_stddevs_train &lt;- attr(train_data, "scaled:scale")
test_data &lt;- scale(test_data, center = col_means_train, scale = col_stddevs_train)

train_data[1, ] # First training sample, normalized

######################################
# OLS
olstrain = data.frame(cbind(train_labels, train_data))
olstest  = data.frame(cbind(test_labels, test_data))
colnames(olstrain) = c('y', 'CRIM', 'ZN', 'INDUS', 'CHAS', 'NOX', 'RM', 'AGE', 
                       'DIS', 'RAD', 'TAX', 'PTRATIO', 'B', 'LSTAT')
colnames(olstest) = c('y', 'CRIM', 'ZN', 'INDUS', 'CHAS', 'NOX', 'RM', 'AGE', 
                      'DIS', 'RAD', 'TAX', 'PTRATIO', 'B', 'LSTAT')
ols = lm(y~.,data=olstrain)
summary(ols)
olspred = predict(ols, newdata=olstest)

library(Metrics)
mae(test_labels, olspred)

##################################
# Keras

#lrate=0.02
lrate=0.1

# OLS-like results are sensitive to choice of LR and EPOCHS (both must match)

epochs = 100
build_model &lt;- function() {
  model &lt;- keras_model_sequential() %&gt;%

    # OLSisch
    layer_dense(units = 1, input_shape = dim(train_data)[2])

    # Simple NN
    #layer_dense(units = 64, activation = "relu", input_shape = dim(train_data)[2]) %&gt;%
    #layer_dense(units = 1)

    model %&gt;% compile(
      loss = "mse",
      optimizer = optimizer_rmsprop(lr=lrate),
      metrics = list("mean_absolute_error")
    )

  model
}

model &lt;- build_model()
#model %&gt;% summary()

# Fit the model and store training stats
history &lt;- model %&gt;% fit(
  train_data,
  train_labels,
  epochs = epochs,
  validation_split = 0.2,
  verbose = 1
)

c(loss, kmae) %&lt;-% (model %&gt;% evaluate(test_data, test_labels, verbose = 0))

print("======================================================")
print(paste0("Learning rate: ", lrate))
print(paste0("Keras: Mean absolute error on test set: $", sprintf("%.2f", kmae*1000)))
print(paste0("Keras mean val_mae:                     $", sprintf("%.2f", mean(history$metrics$val_mean_absolute_error[30:100])*1000)))
print(paste0("O L S: Mean absolute error on test set: $", sprintf("%.2f", mae(test_labels, olspred)*1000)))
print(paste0("Diff.                                 : $", sprintf("%.2f", kmae*1000 - mae(test_labels, olspred)*1000)))
I would say your way of comparing this is rather incorrect. There are two main problem from what you did those are :


R linear model (lm) is solved using deterministic solution whereas
keras model optimized with gradient descent. In the end you expect
the weights and result to be close but not similar. For OLS, it
attains the lowest theoretical training error(by definition)
averaged over all training sample. Your keras model will converge
towards this minima but it might have difficulty converging to
similar minima(due to the nature of GD).
you are not training your keras model with full data. Notice you put
validation split to be 0.2 i.e. you are training with 80 percent of
your data.


Finally, if you want to test whether single layer NN is similar to OLS You can try setting validation split to 0 (or simply don't declare it) and training it with the whole dataset as your batch.