From Data to Predictions: Understanding Linear Regression Through Statistics and Python

Why Statistics Matter in Machine Learning

Machine Learning is often described as "teaching computers to learn from data", but under the hood, most of that learning is really applied statistics. Every time a model makes a prediction, it is quietly asking questions like: what is the typical value in this data? How spread out is it? Do these two things move together?

If you can answer those three questions — the average, the spread, and the relationship — you already have the raw ingredients to build one of the most useful ML models ever invented: Linear Regression. In this article we will walk that path end to end, using one running example: predicting an employee's salary from their years of experience.

• Statistics gives us a language to describe data (mean, variance, covariance).
• It helps us decide which features matter and which are noise.
• It powers the formulas that estimate model coefficients.
• It provides the metrics we use to judge whether a model is any good.

Tip — you do not need a maths degree to follow along. If you can compute an average on paper, you can understand everything that follows. We introduce each formula only after you have felt the intuition first.

Meet the Training Data

A training dataset is simply a table of past examples the model learns from. Each row is one example. Each column is either a feature (an input the model is allowed to see) or the target (the answer we want the model to predict).

For our salary example we will use a tiny, hand-picked dataset of six employees. It is small enough that we can compute every statistic by hand, and big enough to show the intuition. The feature X is Years of Experience. The target Y is Salary in thousands of dollars.

Employee   Years (X)   Salary in $k (Y)
--------   ---------   ----------------
   E1          1              40
   E2          2              50
   E3          3              60
   E4          4              65
   E5          5              75
   E6          6              85

Everyday analogy — think of this as a mentor's notebook. Over the years the mentor has watched juniors grow. Every time someone got a raise, the mentor scribbled down (their years, their new salary). Now a new hire walks in and asks: "Based on your notebook, what should I expect to earn in three years?" Linear Regression is a formal way to answer that question.

import numpy as np
import pandas as pd

# Our tiny salary dataset
data = {
    "YearsExperience": [1, 2, 3, 4, 5, 6],
    "Salary":          [40, 50, 60, 65, 75, 85],  # in $ thousands
}
df = pd.DataFrame(data)
print(df)
print("\nShape:", df.shape)

The output shows six rows and two columns. This is the entire universe the model is allowed to learn from — nothing else. Everything we compute from here on (mean, variance, covariance, slope, intercept) comes from these twelve numbers.

The Mean — Where the Data Sits

The mean (arithmetic average) is the single number that best represents the "centre" of a column. Add up all the values, divide by how many there are, and you have it. It answers the question: if I had to summarise this whole column with one number, what would it be?

For our six employees, the mean of Years is (1+2+3+4+5+6)/6 = 3.5 years, and the mean of Salary is (40+50+60+65+75+85)/6 = 62.5 (that is $62,500). The mean matters for Linear Regression because the best-fit straight line is guaranteed to pass through the point (mean of X, mean of Y). This single fact will pin down our intercept later.

import numpy as np

X = np.array([1, 2, 3, 4, 5, 6])
Y = np.array([40, 50, 60, 65, 75, 85])

mean_x = X.mean()
mean_y = Y.mean()

print(f"Mean of Years  (X-bar) = {mean_x}")
print(f"Mean of Salary (Y-bar) = {mean_y}")

Note — the mean is sensitive to outliers. If one senior with 20 years of experience earning $300k sneaked into our table, the mean salary would jump dramatically even though most employees earn much less. This is why data cleaning matters before modelling.

Variance — How Spread Out the Data Is

Two datasets can share the same mean but look completely different. Imagine two teams both averaging 62.5k in salary: in one team everyone earns close to 62.5k; in the other, half earn 30k and half earn 95k. The mean hides that difference. Variance is the statistic that reveals it — it measures how far, on average, the values are from the mean.

We compute variance in three steps: subtract the mean from every value (that gives the "error" from the centre), square each error (so positives and negatives don't cancel), and then average those squares.

Everyday analogy — variance is like asking "how noisy is this class?" A silent classroom has near-zero variance in decibels. A noisy one has high variance. In Linear Regression, variance of X tells the algorithm how much of a "lever" it has to work with: if X barely changes across the rows, the model cannot learn much from it.

import numpy as np

X = np.array([1, 2, 3, 4, 5, 6], dtype=float)
Y = np.array([40, 50, 60, 65, 75, 85], dtype=float)

mean_x = X.mean()
mean_y = Y.mean()

# Manual (population) variance — sum of squared deviations, divided by n
var_x = ((X - mean_x) ** 2).sum() / len(X)
var_y = ((Y - mean_y) ** 2).sum() / len(Y)

print(f"Variance of X = {var_x:.3f}")
print(f"Variance of Y = {var_y:.3f}")

# NumPy shortcut (population variance uses ddof=0)
print("NumPy check   :", np.var(X), np.var(Y))

You should see Var(X) = 2.917 and Var(Y) ≈ 227.083. Salary varies much more than years — that is expected: salary is in thousands of dollars, years is a small integer. Units always affect variance, which is why in reporting we often prefer the square root of variance (the standard deviation).

Covariance — Do Two Things Move Together?

Mean and variance describe one column at a time. Covariance describes two columns at once: does Y tend to go up when X goes up? If yes, covariance is positive. If Y tends to go down when X goes up, covariance is negative. If there is no pattern, covariance is around zero.

For each row we compute how far X is from its mean, how far Y is from its mean, and multiply. If both are above their means (or both below), the product is positive. If one is above and the other below, it is negative. Averaging those products across all rows gives the covariance.

Everyday analogy — covariance is like watching two dancers. If they lift their arms at the same time, they "co-vary" positively. If one raises while the other lowers, they co-vary negatively. If their movements look random, covariance is near zero. Linear Regression is essentially a formal way of saying: how much does Y dance along with X?

import numpy as np

X = np.array([1, 2, 3, 4, 5, 6], dtype=float)
Y = np.array([40, 50, 60, 65, 75, 85], dtype=float)

mean_x, mean_y = X.mean(), Y.mean()

# Manual (population) covariance
cov_xy = ((X - mean_x) * (Y - mean_y)).sum() / len(X)
print(f"Covariance(X, Y) = {cov_xy:.3f}")

# NumPy check — np.cov returns a 2x2 matrix; use bias=True for population form
print("NumPy check      :", np.cov(X, Y, bias=True)[0, 1])

You should see Cov(X, Y) ≈ 25.417 — a strong positive number that confirms our intuition: as years of experience go up, salary goes up. That single number is about to become the engine of our model.

From Statistics to Coefficients

Linear Regression tries to draw the straight line that best passes through our scatter plot of (X, Y) points. Every straight line can be written as Y = b0 + b1 * X, where b1 is the slope (how much Y changes when X goes up by 1) and b0 is the intercept (the value of Y when X is 0). Our job is to find the two numbers b0 and b1 that make the line fit the data as closely as possible.

It turns out that the best line — the one that minimises the total squared error between the points and the line — has beautifully simple formulas. The slope is the covariance of X and Y divided by the variance of X. The intercept is fixed by the requirement that the line must pass through the point of means.

Plug in our numbers: b1 = 25.417 / 2.917 ≈ 8.714. So every extra year of experience is worth about $8,714 in expected salary. Then b0 = 62.5 - 8.714 * 3.5 ≈ 31.999, meaning a brand-new employee with zero years of experience is estimated to earn around $32,000. Notice how mean, variance, and covariance all showed up in a single, elegant recipe.

import numpy as np

X = np.array([1, 2, 3, 4, 5, 6], dtype=float)
Y = np.array([40, 50, 60, 65, 75, 85], dtype=float)

mean_x, mean_y = X.mean(), Y.mean()

cov_xy = ((X - mean_x) * (Y - mean_y)).sum() / len(X)
var_x  = ((X - mean_x) ** 2).sum() / len(X)

b1 = cov_xy / var_x           # slope
b0 = mean_y - b1 * mean_x     # intercept

print(f"Slope     (b1) = {b1:.4f}")
print(f"Intercept (b0) = {b0:.4f}")
print(f"Model: Salary = {b0:.2f} + {b1:.2f} * YearsExperience")

Common mistake — students sometimes divide the slope by variance of Y instead of variance of X. Always divide by the variance of the feature (the input), not the target. Intuitively, we are asking: "per unit of X, how much does Y move?", so the denominator must be in the units of X.

Making Predictions With the Trained Line

Once we have b0 and b1, prediction is arithmetic. For any new value of X we simply compute Y_hat = b0 + b1 * X. The hat symbol reminds us this is the model's estimate, not the true (unknown) salary.

b0 = 31.9999
b1 = 8.7143

def predict(years):
    return b0 + b1 * years

for yrs in [0, 3, 7, 10]:
    print(f"{yrs:2d} years  ->  predicted salary = $ {predict(yrs):.2f}k")

The model estimates about $58k for 3 years, $93k for 7 years, and $119k for 10 years. Note the graceful extrapolation — but be careful. Predicting far outside the range of the training data (say, 40 years of experience) is unreliable because we have no evidence the straight-line pattern still holds that far out.

Building the Same Model With Scikit-learn

In production you rarely compute coefficients by hand — you let a well-tested library do it. Scikit-learn's LinearRegression class implements exactly the same math we did above (plus a lot of numerical safety). If our manual work is correct, the two answers should match to several decimal places.

import numpy as np
from sklearn.linear_model import LinearRegression

X = np.array([[1], [2], [3], [4], [5], [6]], dtype=float)  # 2D shape
Y = np.array([40, 50, 60, 65, 75, 85], dtype=float)

model = LinearRegression()
model.fit(X, Y)

print(f"Scikit-learn intercept (b0) = {model.intercept_:.4f}")
print(f"Scikit-learn slope     (b1) = {model.coef_[0]:.4f}")

# Predict for the same inputs used in the manual example
for yrs in [0, 3, 7, 10]:
    pred = model.predict([[yrs]])[0]
    print(f"{yrs:2d} years  ->  predicted salary = $ {pred:.2f}k")

You should see b0 ≈ 31.9999 and b1 ≈ 8.7143 — identical to our hand calculation. Note that scikit-learn expects X to be a 2D array (rows × features) even when there is only one feature. This tiny detail trips up almost every beginner at least once.

Judging the Model — MAE, MSE, RMSE, and R²

How do we know if our line is any good? We compare its predictions against the actual salaries in the training data (or, better, a held-out test set). The gap between predicted and actual is called the residual. Different metrics summarise those residuals in different ways — each with its own personality.

• MAE (Mean Absolute Error) — average of |actual - predicted|. Same units as Y. Easy to explain: "on average we are off by X dollars".
• MSE (Mean Squared Error) — average of (actual - predicted)². Punishes big misses more than small ones. Units are Y-squared, so harder to interpret.
• RMSE (Root Mean Squared Error) — square root of MSE. Brings the units back to Y. The most reported metric in industry.
• R² Score — the fraction of the variance in Y that the model successfully explains. 1.0 is perfect, 0.0 means the model does no better than predicting the mean.

import numpy as np
from sklearn.linear_model import LinearRegression
from sklearn.metrics import mean_absolute_error, mean_squared_error, r2_score

X = np.array([[1], [2], [3], [4], [5], [6]], dtype=float)
Y = np.array([40, 50, 60, 65, 75, 85], dtype=float)

model = LinearRegression().fit(X, Y)
Y_pred = model.predict(X)

mae  = mean_absolute_error(Y, Y_pred)
mse  = mean_squared_error(Y, Y_pred)
rmse = np.sqrt(mse)
r2   = r2_score(Y, Y_pred)

print(f"MAE  = {mae:.3f} (thousand $)")
print(f"MSE  = {mse:.3f}")
print(f"RMSE = {rmse:.3f} (thousand $)")
print(f"R^2  = {r2:.4f}")

For this dataset R² is about 0.987, meaning the straight line explains ~98.7% of the variation in salary — an excellent fit for such a simple model. Interpreting in business terms: "On average our salary estimate is within ~$1.5k of the truth (MAE), and knowing an employee's years of experience alone explains almost all of the salary variation across the team."

Warning — a very high R² on training data does not automatically mean a great model. Always evaluate on data the model has never seen (a test set), otherwise you may be measuring memorisation instead of learning.

Seeing the Line — Scatter Plot and Fit

A picture will make everything above click. Let's plot the six real employees as dots and overlay the regression line. Ideally the dots should sit close to the line, with the line passing through the centre of the cloud.

import numpy as np
import matplotlib.pyplot as plt
from sklearn.linear_model import LinearRegression

X = np.array([[1], [2], [3], [4], [5], [6]], dtype=float)
Y = np.array([40, 50, 60, 65, 75, 85], dtype=float)

model = LinearRegression().fit(X, Y)

# A smooth line from 0 to 8 years for plotting
xs = np.linspace(0, 8, 100).reshape(-1, 1)
ys = model.predict(xs)

plt.figure(figsize=(7, 4.5))
plt.scatter(X, Y, s=70, color="#2563eb", label="Actual salaries")
plt.plot(xs, ys, color="#e11d48", linewidth=2, label="Regression line")
plt.scatter([X.mean()], [Y.mean()], s=140, marker="X",
            color="black", label="(mean X, mean Y)")
plt.xlabel("Years of Experience")
plt.ylabel("Salary ($k)")
plt.title("Salary vs Years of Experience — Linear Regression Fit")
plt.legend()
plt.grid(alpha=0.3)
plt.show()

Notice the black X marker in the middle — that is the point (mean X, mean Y) = (3.5, 62.5). The regression line must pass through it. This is the visual confirmation of the formula b0 = mean(Y) - b1 * mean(X).

Applying It to a Larger Real-World Dataset

Six rows is enough to build intuition but too small for a serious model. Let's now scale up to a realistic "Years of Experience vs Salary" dataset with 30 employees — the same shape of problem, but with the natural noise you would see in a real HR database.

import numpy as np
import pandas as pd
from sklearn.linear_model import LinearRegression
from sklearn.model_selection import train_test_split
from sklearn.metrics import mean_absolute_error, mean_squared_error, r2_score

# A realistic (synthetic) salary dataset — 30 employees
rng = np.random.default_rng(42)
years  = np.round(rng.uniform(0.5, 12, size=30), 1)
noise  = rng.normal(0, 4, size=30)          # real-world noise
salary = 30 + 8.5 * years + noise           # true rule + noise

df = pd.DataFrame({"YearsExperience": years, "Salary": np.round(salary, 2)})
print(df.head(10))

# 1) Split into train / test so we evaluate on unseen data
X = df[["YearsExperience"]].values
y = df["Salary"].values
X_train, X_test, y_train, y_test = train_test_split(
    X, y, test_size=0.25, random_state=1
)

# 2) Fit the model on training data only
model = LinearRegression().fit(X_train, y_train)
print(f"\nLearned model:  Salary = {model.intercept_:.2f} + "
      f"{model.coef_[0]:.2f} * YearsExperience")

# 3) Evaluate on the held-out test set
y_pred = model.predict(X_test)
print(f"MAE  = {mean_absolute_error(y_test, y_pred):.2f} k$")
print(f"RMSE = {np.sqrt(mean_squared_error(y_test, y_pred)):.2f} k$")
print(f"R^2  = {r2_score(y_test, y_pred):.4f}")

# 4) Predict salary for three hypothetical new hires
for yrs in [1.5, 5, 9]:
    print(f"Expected salary at {yrs} yrs  = $ {model.predict([[yrs]])[0]:.2f}k")

The learned intercept and slope should be very close to the true values (30 and 8.5) that we baked into the synthetic data — Linear Regression has "discovered the rule" from noisy examples. In a real project you would follow the exact same recipe on real HR data, replacing the synthetic generator with pd.read_csv("salaries.csv").

Business interpretation — if HR asks "what should we offer a candidate with 5 years of experience?", the answer becomes data-driven instead of a gut feel. If HR then asks "how confident are you?", the RMSE gives a plain-English answer: "about ± that many thousand dollars on typical hires".

Metric Cheat Sheet

Metric   Same units as Y?   Sensitive to outliers?   Best when...
------   ----------------   ----------------------   ---------------------------
MAE      Yes                Low                       You want a plain-English error
MSE      No (squared)       Very high                 Used inside training/optimisation
RMSE     Yes                High                      You report to stakeholders
R^2      Unitless (0..1)    Medium                    You compare models on same target

Common Interview Questions

Beginner Level (10 Questions)

Q1. What is Linear Regression?

Linear Regression is a supervised ML algorithm that models the relationship between one or more input features X and a continuous target Y by fitting a straight line (or hyperplane in higher dimensions) of the form Y = b0 + b1*X. The goal is to learn the coefficients b0 and b1 that minimise the squared error between the predicted and the actual Y values.

Q2. What is the difference between regression and classification?

Both are supervised learning, but regression predicts a continuous number (salary, price, temperature) while classification predicts a discrete label (spam / not-spam, churn / retain). The choice of evaluation metric also differs — RMSE and R² for regression, accuracy and F1 for classification.

Q3. Why do we compute the mean of a variable?

The mean gives us a single, representative value for the column. In Linear Regression it also has a special role: the best-fit line is always guaranteed to pass through the point (mean X, mean Y). This pins the intercept once the slope is known.

Q4. What is variance and why do we square the deviations?

Variance measures how spread out the values are from their mean. We square the deviations so that positive and negative errors don't cancel and so that large deviations are penalised more than small ones. The result is always non-negative.

Q5. What does covariance tell us?

Covariance measures whether two variables move together. Positive covariance means both tend to rise together (like experience and salary). Negative covariance means one falls when the other rises. A covariance near zero suggests no linear relationship.

Q6. How is the slope in Linear Regression estimated?

The slope b1 equals the covariance of X and Y divided by the variance of X: b1 = Cov(X, Y) / Var(X). Intuitively, you are asking "per unit change in X, how much does Y change on average?"

Q7. How is the intercept estimated?

Once the slope is known, the intercept follows from the constraint that the line must pass through the point of means: b0 = mean(Y) - b1 * mean(X). It represents the predicted value of Y when X = 0.

Q8. Why do we need a training and a test set?

The training set teaches the model, the test set independently measures how well the model performs on data it has never seen. Without a test set you can't tell whether the model has learned the underlying pattern or simply memorised the training rows (overfitting).

Q9. What is R² and what does R² = 0.87 mean?

R² (coefficient of determination) is the fraction of the variance of Y that the model explains. R² = 0.87 means about 87% of the variation in the target is explained by the model; the remaining 13% is noise or missing features. R² = 1 is a perfect fit, R² = 0 is no better than predicting the mean, and negative R² means worse than the mean.

Q10. When should we NOT use Linear Regression?

Avoid it when the relationship is clearly non-linear (e.g. exponential growth), when there are strong interactions between features, when the target is categorical, or when the data has heavy outliers that violate the model's assumptions. In those cases consider polynomial regression, tree-based models, or classification algorithms.

Intermediate Level (5 Questions)

Q11. What are the assumptions of Linear Regression?

Linearity (X and Y are linearly related), independence of observations, homoscedasticity (constant variance of residuals), normality of residuals, and low multicollinearity between features. Violating these does not automatically break the model but usually degrades the reliability of its coefficient estimates and p-values.

Q12. What is the cost function that Linear Regression minimises?

The Ordinary Least Squares (OLS) cost function is the mean (or sum) of squared residuals: J(b0, b1) = (1/n) * Σ (y_i - (b0 + b1 * x_i))². The closed-form formulas for b0 and b1 come from setting the partial derivatives of J to zero and solving.

Q13. Difference between MSE and RMSE — and which do we report?

MSE is the average squared error; RMSE is its square root. RMSE has the same units as Y, so it is much easier to explain to stakeholders ("we are typically off by $3k"). MSE is convenient inside optimisation because it is smooth and differentiable, but RMSE is preferred in reports.

Q14. What is multicollinearity and why is it a problem?

Multicollinearity happens when two or more input features are highly correlated with each other (e.g. "years of experience" and "age"). The model can still make good predictions overall, but the individual coefficients become unstable and hard to interpret. Detect it with the Variance Inflation Factor (VIF) and treat it by dropping or combining redundant features.

Q15. Why can R² be misleading, and what alternatives exist?

R² always increases (or stays the same) when you add more features, even useless ones, so it rewards complexity. Adjusted R² penalises unnecessary features. On new data, you should always accompany R² with RMSE or MAE on a test set and, when possible, cross-validated scores.

Hands-on Practice

• Exercise 1 — Extend the six-row salary dataset to ten rows. Recompute mean, variance, and covariance by hand, then confirm with NumPy.
• Exercise 2 — Change one salary from 85 to 200 (a big outlier). Refit the model and observe how much the slope changes. What does this tell you about outliers?
• Exercise 3 — Manually calculate b0 and b1 for the dataset X = [2, 4, 6, 8], Y = [3, 7, 5, 10]. Then verify with scikit-learn.
• Exercise 4 — Download a public "Salary_Data.csv" (Years vs Salary) from Kaggle. Fit a LinearRegression, split 80/20, and report MAE, RMSE, and R² on the test set.
• Exercise 5 — Predict house prices instead. Pick a small housing dataset with a single feature such as area (sq ft). Reuse the same code with only the column names changed.
• Exercise 6 — Plot the residuals (actual - predicted) vs X. If you see a pattern (a curve or a funnel shape), it is a hint that a straight line is not the best model.

import numpy as np
from sklearn.linear_model import LinearRegression
from sklearn.metrics import mean_absolute_error, r2_score

# Change these two lines and rerun to answer Exercise 3
X = np.array([[2], [4], [6], [8]], dtype=float)
Y = np.array([3, 7, 5, 10], dtype=float)

model = LinearRegression().fit(X, Y)
print(f"b0 = {model.intercept_:.4f}")
print(f"b1 = {model.coef_[0]:.4f}")
print(f"R^2 = {r2_score(Y, model.predict(X)):.4f}")
print(f"MAE = {mean_absolute_error(Y, model.predict(X)):.4f}")

Bringing It All Together

We started from a simple question — how do we predict a salary from years of experience — and let statistics guide us all the way to a working model. The mean told us where the data sits. The variance told us how spread out it is. The covariance told us whether X and Y dance together. Combined into two little formulas, they gave us the slope and intercept of the best-fit line, and that line gave us predictions we can actually use in business.

Every more advanced ML model — from multiple linear regression to gradient-boosted trees to deep neural networks — builds on this same statistical foundation. Once you truly feel mean, variance, and covariance in your bones, the rest of ML becomes far less mysterious.

What's Next?

• Multiple Linear Regression — add more features (age, education, city) and see how the formulas generalise to matrices.
• Ordinary Least Squares (OLS) — learn the matrix form b = (XᵀX)⁻¹Xᵀy that scikit-learn actually uses under the hood.
• Gradient Descent — an iterative alternative to the closed-form solution, essential for training large models where matrix inversion is infeasible.
• Polynomial Regression — model curved relationships by adding X², X³ features while still using a linear model.
• Regularization (Ridge, Lasso, Elastic Net) — control overfitting and perform automatic feature selection.
• Logistic Regression — the natural next step when the target becomes a category instead of a number.

Tip — do not rush to complex models. Master Linear Regression on three or four real datasets first. A senior ML engineer who understands why Cov(X,Y) / Var(X) is the slope will out-diagnose a beginner who can only call model.fit().