After having a bit too much tea I decided to write a silly way of estimating the parameters of a single-feature linear regression. For starters, I avoided using techniques like gradient descent, least-squares, or a normal equation. I also didn't want to make the regression any more sophisticated by using feature scaling or regularization. Indeed, vectorization was intentionally avoided as well.
Despite my efforts to thwart validity, I think what I conjured up offers a suboptimal and inefficient, but not entirely ineffective, way to parametrize a linear model. I was able to fit it to some randomly generated data almost as well as least-squares, although I did find that it was outlier sensitive (which is expected from using averages).
Essentially, I used the two-point slope calculation for each data point with each subsequent datapoint, and took the arithmetic mean to estimate the slope. The error in the population slope was estimated with the standard deviation of the two-point slopes. Using each calculated slope with each data point, I calculate the intercept. The arithmetic mean of the calculated intercepts is used to estimate the population intercept, and the error was again estimated with the standard deviation.
As I see it, this is basically descriptive statistics done on multiple two-point estimates of parameters, but I want to know what others think.
What does this method of parametrization do invalidly, and validly? I would like to know what this method has been called if it has a name. What assumptions does this method make?
I've supplied the implementation for others to better understand the method, and to try it out if they want.
import numpy as np import matplotlib.pyplot as plt def regression(x, y, plot=False): slopes =  for i in range(len(y)): for j in range(len(y)): if (j>i): slopes.append((y[j] - y[i]) / (x[j] - x[i])) slopes_mean = np.mean(slopes) slopes_std = np.std(slopes) intercepts =  for slope in slopes: for i in range(len(y)): intercepts.append(y[i] - slope * x[i]) intercepts_mean = np.mean(intercepts) intercepts_std = np.std(intercepts) residuals =  for slope in slopes: for intercept in intercepts: for i in range(len(y)): residuals.append((intercept + slope * x[i] - y[i])**2) residual_error = sum(residuals) if plot == True: plt.scatter(x,y,color='k') plt.plot([i for i in np.arange(min(x),max(x) + 1)],[i*np.mean(slopes) + np.mean(intercepts) for i in np.arange(min(x),max(x) + 1)], color='r') plt.xlabel('x') plt.ylabel('y') plt.show() return slopes_mean, slopes_std, intercepts_mean, intercepts_std, residual_error
The values of y are predicted exactly by the model when the sample points make a straight line.
The estimated parameters by this method were not much different than those parametrized by scipy.stats.linregress().
The estimated parameters were not far off from the true parameters of a simulated dataset with Gaussian noise added.