PYTHON / MACHINE LEARNING WITH PYTHON
Linear regression
Fit, read and sanity-check a least squares linear regression in scikit-learn, and say what each coefficient and the R^2 score actually mean.
What you will learn
- Shape X as a 2D (n_samples, n_features) array before calling fit
- Read intercept_ and coef_ as 'predicted change in y per unit of that feature'
- Recognise that 'linear' constrains the coefficients, not the features
- Spot under-fitting from curved residuals instead of trusting a single R^2
Understanding Linear regression
Linear regression assumes the target is a weighted sum of the features plus a constant: y = b + w1*x1 + w2*x2 + ... . Fitting means choosing b and the w values that make the sum of squared vertical gaps between the predictions and the observed y as small as possible. Squared errors are used because the resulting objective is smooth and convex in the weights, so there is a single best answer that can be computed directly by linear algebra rather than searched for. That is why LinearRegression has no learning rate, no number of iterations and no random seed: two runs on the same data give bit-for-bit the same weights.
Each coefficient is a partial slope. A coefficient of 3 on x1 means that if x1 rises by one of its own units while every other feature in the model stays fixed, the prediction rises by 3. Two consequences follow immediately: the coefficient's size depends on the unit you measured that feature in, so a coefficient of 0.001 can matter more than one of 500; and the value depends on which other features are in the model, because 'holding the others fixed' changes when you add or drop a column.
The word linear refers to the parameters, not to the shape of the curve. Nothing stops you from handing the model x and x**2 as two separate columns; the model is still a weighted sum, so it still solves in closed form, but the fitted function through the original x is a parabola. This is the mental model to keep: linear regression fits a flat surface in whatever feature space you give it, and your job is to give it a space where a flat surface is a reasonable description. When it isn't, the model fails quietly, with residuals that curve instead of scattering randomly.
R^2 reports the fraction of the target's variance the model accounts for, so 0 means 'no better than always predicting the mean' and 1 means every residual is zero. On four hand-made points that fit an equation exactly you will get 1.0, which says the algebra worked, not that the model would predict anything new.
import numpy as np
from sklearn.linear_model import LinearRegression
# y is built to be exactly 5 + 3*x1 + 2*x2, so least squares must recover it
X = np.array([[1.0, 0.0],
[0.0, 1.0],
[2.0, 1.0],
[3.0, 4.0]])
y = 5 + 3 * X[:, 0] + 2 * X[:, 1]
model = LinearRegression().fit(X, y)
print("intercept:", round(float(model.intercept_), 6))
print("coefficients:", [round(float(c), 6) for c in model.coef_])
print("R^2 on the training data:", round(float(model.score(X, y)), 6))
print("predict x1=10, x2=10:", round(float(model.predict([[10.0, 10.0]])[0]), 6))Linear regression solves in closed form for the weights that minimise squared vertical error, which makes each weight the predicted change in y per unit of that feature with the others held fixed.
Worked examples
Linear in the weights, curved in x
Adding x**2 as a second column lets the same linear model fit a quadratic exactly.
import numpy as np
from sklearn.linear_model import LinearRegression
x = np.array([0.0, 1.0, 2.0, 3.0, 4.0])
y = 2.0 * x**2 - x + 1.0 # exactly quadratic, no noise
X1 = x.reshape(-1, 1) # just x
X2 = np.column_stack([x, x**2]) # x and x squared
straight = LinearRegression().fit(X1, y)
curved = LinearRegression().fit(X2, y)
print("R^2 with x only:", round(float(straight.score(X1, y)), 4))
print("R^2 with x and x**2:", round(float(curved.score(X2, y)), 4))
print("curved intercept:", round(float(curved.intercept_), 4))
print("curved coefficients:", [round(float(c), 4) for c in curved.coef_])Example explained
Line 1x.reshape(-1, 1) turns 5 values into 5 rows of 1 feature, which is the shape fit expects.
Line 2np.column_stack([x, x**2]) builds a 5x2 design matrix; the model treats x**2 as an unrelated feature.
Line 3R^2 of 0.8974 looks respectable even though a straight line cannot possibly fit a parabola, which is why a single score is weak evidence.
Line 4The recovered weights -1.0 and 2.0 with intercept 1.0 are exactly the coefficients of 2*x**2 - x + 1.
What fit is actually computing
Solving the least squares problem yourself with np.linalg.lstsq reproduces sklearn's intercept and slope.
import numpy as np
from sklearn.linear_model import LinearRegression
X = np.array([[1.0], [2.0], [3.0], [4.0]])
y = np.array([2.0, 4.1, 6.2, 7.9])
# a column of ones makes the intercept just one more coefficient
A = np.hstack([np.ones((4, 1)), X])
beta, *_ = np.linalg.lstsq(A, y, rcond=None)
sk = LinearRegression().fit(X, y)
print("lstsq intercept, slope:", [round(float(b), 6) for b in beta])
print("sklearn intercept, slope:", [round(float(sk.intercept_), 6),
round(float(sk.coef_[0]), 6)])Example explained
Line 1np.hstack adds the constant column, so the intercept is fitted the same way as any weight.
Line 2np.linalg.lstsq returns the weight vector that minimises the squared residuals, plus diagnostics that are discarded here.
Line 3The two printed lines agree because LinearRegression is a wrapper over exactly this least squares solve, not an iterative optimiser.
Line 4The data is noisy, so neither result is a perfect fit; both find the same single best compromise.
Important notes
Because errors are squared, one far-off point pulls the fitted line much harder than several small errors combined, so inspect outliers before blaming the model.
If two columns are exact duplicates or one is a sum of others, sklearn's SVD-based solver does not raise an error; it returns one arbitrary split of the shared effect, so those individual coefficients are not interpretable.
Common mistakes
Passing a 1D array such as np.array([1, 2, 3]) as X: fit raises ValueError telling you to reshape with .reshape(-1, 1), because one row per sample cannot be inferred from a flat list.
Ranking features by the size of coef_ when the columns are in different units, for example metres against millimetres, which makes an important feature look negligible purely because its unit is large.
Treating a high training R^2 as proof the relationship is linear: a straight line through clearly curved data still scores well, and the only signal is the residual pattern.
Try it yourself
Change, predict, then run
Fit LinearRegression to X = [[1], [2], [3], [4], [5]] and y = [3.1, 5.0, 6.9, 9.2, 11.0], then print the intercept and coefficient and confirm they are close to 1 and 2. Predict the value at x = 10 and check it lands near 21.
Open the Python workspaceCheck your understanding
A model predicts house price from area measured in square metres and gets a coefficient of 3000. You refit the identical model with area expressed in square centimetres instead. What happens?
- The coefficient stays 3000 and every prediction becomes 10000 times larger
- The coefficient becomes 0.3 and the predictions and R^2 are unchanged
- The coefficient becomes 30000000 and R^2 improves because the feature carries more detail
- The coefficient stays 3000 and R^2 drops because the feature values are now very large
Show answer
A coefficient is priced per one unit of the feature, and one square metre is 10000 square centimetres, so the weight shrinks by the same factor and coefficient times feature value is identical. The fitted function, the residuals and therefore R^2 do not move at all. The last option is tempting because scaling genuinely matters for gradient-based or regularised fitting, but plain least squares is solved exactly and a change of unit only relabels the same line.