Compute all pairwise Pearson r values for three columns and assemble the 3×3 symmetric correlation matrix (diagonal = 1). Use sum() for column means, then one joint loop accumulating six sums (three cross-products and three squared-deviation sums). Three distinct off-diagonal r values reveal that each pair has a different relationship: strong positive (a,b), moderate negative (a,c), weak negative (b,c). With numpy, np.corrcoef([a, b, c]).

By hand

With the library

np.corrcoef([a, b, c]) returns the full 3×3 matrix in one call. Diagonals are always 1.0; the three unique off-diagonal values are m[0,1], m[0,2], m[1,2]. Note the three relationships are visibly distinct — not all equal.

naive.py
import math
a = [1, 2, 3, 4]
b = [1, 2, 4, 3]
c = [4, 1, 3, 2]
n = len(a)
ma = sum(a) / n
mb = sum(b) / n
mc = sum(c) / n
ss_ab = 0.0
ss_ac = 0.0
ss_bc = 0.0
ss_a = 0.0
ss_b = 0.0
ss_c = 0.0
for i in range(n):
    ss_ab = ss_ab + (a[i] - ma) * (b[i] - mb)
    ss_ac = ss_ac + (a[i] - ma) * (c[i] - mc)
    ss_bc = ss_bc + (b[i] - mb) * (c[i] - mc)
    ss_a = ss_a + (a[i] - ma) ** 2
    ss_b = ss_b + (b[i] - mb) ** 2
    ss_c = ss_c + (c[i] - mc) ** 2
r_ab = ss_ab / math.sqrt(ss_a * ss_b)
r_ac = ss_ac / math.sqrt(ss_a * ss_c)
r_bc = ss_bc / math.sqrt(ss_b * ss_c)
print('RESULT:', (round(r_ab, 4), round(r_ac, 4), round(r_bc, 4)))
library.py
import numpy as np
from dalib.display import set_display
set_display()

a = [1, 2, 3, 4]
b = [1, 2, 4, 3]
c = [4, 1, 3, 2]
m = np.corrcoef([a, b, c])
print('r_ab:', round(float(m[0, 1]), 4))
print('r_ac:', round(float(m[0, 2]), 4))
print('r_bc:', round(float(m[1, 2]), 4))
print('RESULT:', (round(float(m[0, 1]), 4), round(float(m[0, 2]), 4), round(float(m[1, 2]), 4)))
r_ab: 0.8
r_ac: -0.4
r_bc: -0.2
RESULT: (0.8, -0.4, -0.2)

Honesty

This lesson shows the computation exactly, on a tiny pinned sample. The pairwise correlations are correct and reproducible, but with this few rows they demonstrate the mechanism, not evidence — a correlation matrix on a handful of observations is not a valid statistical finding. Real inference needs an adequate sample size and assumption checks (independence, linearity, bivariate normality); read each entry as "how the coefficient is computed," not as a conclusion about the population.

Implementation notes

  • sum(col)/n treats sum() as a primitive (mean computation was shown in pearson-correlation; here the focus is the cross-product accumulation).
  • The six-accumulator joint loop computes all three pairwise r values in a single O(n) pass, avoiding three separate paired loops.
  • Interpreting the matrix: (a,b)=0.8 — b tracks a closely (they differ only in positions 2 and 3); (a,c)=−0.4 — c partially inverts a; (b,c)=−0.2 — weak negative, little linear structure.
  • For k columns the number of unique off-diagonal pairs is k(k−1)/2 (here k=3 gives 3 pairs). np.corrcoef handles arbitrary k in one call.
  • Cross-reference: pearson-correlation (this chapter) for the single-pair formula that each matrix cell applies.