Measure how two variables move together: sample covariance = sum((xi−mx)× (yi−my)) / (n−1). Three loops: one for mean_x, one for mean_y, one for the sum of products. Divide by n−1 (sample covariance, matching the variance convention). With numpy, np.cov(x, y)[0, 1] defaults to ddof=1 and gives the same value. Data: n=5 pairs, mx=3, my=4, cov_sum=8, cov=2.0.

By hand

Loop 1: accumulate total_x → mx=3.0. Loop 2: accumulate total_y → my=4.0. Loop 3: for each index i, accumulate (x[i]−mx)×(y[i]−my) into cov_sum. Deviations: x=[-2,-1,0,1,2], y=[-2,0,0,0,2]; products=[4,0,0,0,4], cov_sum=8. Divide by n−1=4: cov=2.0.

naive.py
Replay: real traced execution (multi-file project)
x = [1, 2, 3, 4, 5]
y = [2, 4, 4, 4, 6]
n = len(x)
total_x = 0.0
for v in x:
    total_x = total_x + v
mx = total_x / n
total_y = 0.0
for v in y:
    total_y = total_y + v
my = total_y / n
cov_sum = 0.0
for i in range(n):
    cov_sum = cov_sum + (x[i] - mx) * (y[i] - my)
cov = cov_sum / (n - 1)
print('RESULT:', round(cov, 10))
  1. x ← [1, 2, 3, 4, 5]

    1x = [1, 2, 3, 4, 5]2y = [2, 4, 4, 4, 6]
    values this step[1, 2, 3, 4, 5]x
  2. y ← [2, 4, 4, 4, 6]

    1x = [1, 2, 3, 4, 5]2y = [2, 4, 4, 4, 6]3n = len(x)
    values this step[2, 4, 4, 4, 6]y
  3. n ← 5

    2y = [2, 4, 4, 4, 6]3n = len(x)4total_x = 0.0
    values this step5n
  4. total_x ← 0.0

    3n = len(x)4total_x = 0.05for v in x:
    values this step0.0total_x
  5. v ← 1, total_x ← 1.0

    pass 1 of 5
    4total_x = 0.05for v in x:6    total_x = total_x + v7mx = total_x / n
    values this step1v0.0 1.0total_x
    All 5 passes — pass 1 is the card above
    passvtotal_x
    110.0 1.0
    21 21.0 3.0
    32 33.0 6.0
    43 46.0 10.0
    54 510.0 15.0
  6. for v in x:

    4total_x = 0.05for v in x:6    total_x = total_x + v
  7. mx ← 3.0

    6    total_x = total_x + v7mx = total_x / n8total_y = 0.0
    values this step3.0mx
  8. total_y ← 0.0

    7mx = total_x / n8total_y = 0.09for v in y:
    values this step0.0total_y
  9. v ← 2, total_y ← 2.0

    pass 1 of 5
    8total_y = 0.09for v in y:10    total_y = total_y + v11my = total_y / n
    values this step5 2v0.0 2.0total_y
    All 5 passes — pass 1 is the card above
    passvtotal_y
    15 20.0 2.0
    22 42.0 6.0
    36.0 10.0
    410.0 14.0
    54 614.0 20.0
  10. for v in y:

    8total_y = 0.09for v in y:10    total_y = total_y + v
  11. my ← 4.0

    10    total_y = total_y + v11my = total_y / n12cov_sum = 0.0
    values this step4.0my
  12. cov_sum ← 0.0

    11my = total_y / n12cov_sum = 0.013for i in range(n):
    values this step0.0cov_sum
  13. i ← 0, cov_sum ← 4.0

    pass 1 of 5
    12cov_sum = 0.013for i in range(n):14    cov_sum = cov_sum + (x[i] - mx) * (y[i] - my)15cov = cov_sum / (n - 1)
    values this step0i0.0 4.0cov_sum
    All 5 passes — pass 1 is the card above
    passicov_sum
    100.0 4.0
    20 1
    31 2
    42 3
    53 44.0 8.0
  14. for i in range(n):

    12cov_sum = 0.013for i in range(n):14    cov_sum = cov_sum + (x[i] - mx) * (y[i] - my)
  15. cov ← 2.0

    14    cov_sum = cov_sum + (x[i] - mx) * (y[i] - my)15cov = cov_sum / (n - 1)16print('RESULT:', round(cov, 10))
    values this step2.0cov
  16. stdout ← RESULT: 2.0

    15cov = cov_sum / (n - 1)16print('RESULT:', round(cov, 10))
    values this stepRESULT: 2.0stdout

With the library

np.cov(x, y) returns a 2×2 covariance matrix; [0, 1] extracts the off-diagonal element (covariance of x and y). Default ddof=1 — same denominator as the naive loop.

library.py
import numpy as np
from dalib.display import set_display
set_display()

x = [1, 2, 3, 4, 5]
y = [2, 4, 4, 4, 6]
cov = np.cov(x, y)[0, 1]
print('cov (ddof=1):', round(float(cov), 10))
print('RESULT:', round(float(cov), 10))
cov (ddof=1): 2.0
RESULT: 2.0

Honesty

This lesson shows the computation of covariance exactly, on a tiny pinned sample. The value is correct and reproducible, but on so few points it demonstrates the mechanism, not a population conclusion — the sign and magnitude here are not a valid statistical finding. Real inference needs an adequate sample size and assumption checks; read it as "how the statistic is computed," not as evidence about how the two variables move together in general.

Implementation notes

  • Sign indicates direction: positive means x and y tend to rise together; negative means one rises as the other falls; zero means no linear relationship.
  • Covariance is scale-dependent: doubling all y values doubles the covariance. This makes comparing covariances across datasets unintuitive — normalizing by the standard deviations gives the unitless Pearson r.
  • Sample (n−1) vs population (n): use n−1 to get an unbiased estimate of the population covariance from a sample, consistent with sample-variance.
  • Cross-reference: pearson-correlation (this chapter) for the normalized version; sample-variance (ch02) for the n−1 convention.