Margin-Based Classifier

The previous two families based their decision on the whole training set with equal weight: in the neighbor procedure each of the 31 neighbors cast a vote, in the counting family each row grew a cell by one unit.

This lesson takes a family that does the opposite. The margin-based classifier places the boundary in the middle of the gap between the two classes and concerns itself only with the examples that fall into that gap; the width of the gap is called the margin. The family also carries a construct that can produce a curved boundary without opening columns: the kernel transform. The lesson counts all three: the relation between the margin width and the held-out numbers, the number of examples that carry the decision, and what the kernel gains against explicit columns.

• NM17. The set, the split and the baseline do not change: a fictional table of 1,260 subscribers, seed 20260218, split 756/252/252, baseline 0.7579 on the test set and 0.8135 on validation. The feature set is the twelve columns that come from the Data Preparation and Feature Engineering course; the scale is learned from the training set alone.
• NM18. The labels are written as $+1$ and $-1$. The model minimizes the sum of the margin loss and the squared penalty on the weights; the penalty coefficient is the tuning variable of this lesson and the margin width is $2 / \lVert w \rVert$.
• NM19. A support vector is a training row that stays inside the margin or on the wrong side of it; those outside the margin contribute nothing to the loss.
• NM20. The solver is a subgradient descent whose step shrinks as $1/(\lambda t)$ and runs 300 rounds in the primal form and 40 rounds in the dual form; the numbers printed are the numbers of this fixed budget.
• NM21. A kernel is a function that gives the inner product in an expanded space without building that space. The polynomial kernel is $(x \cdot z + 1)^2/10$ and the radial kernel is $e^{-0.05 \lVert x - z \rVert^2}$; the decimal coefficients are fixed so that the same penalty range can be used with all three kernels. The dual form keeps as many kernel values as the square of the number of training rows.
• NM22. The choice is read from the validation set; in a tie the wider margin is chosen.

The Width of the Margin

The setup generates the set, writes the model and sweeps the penalty coefficient over eight values.

# margin.py — MODEL: the FICTIONAL subscriber table of the Introduction to Machine Learning course is rebuilt with the same seed;
# the feature set comes from the Data Preparation and Feature Engineering course and is fixed throughout the course.
import math
import statistics

SEED, RAW, M32 = 20260218, 1400, 0xFFFFFFFF
REGION = [("north", 0.28, 21, 0.00), ("south", 0.22, 17, -0.10),
          ("east", 0.18, 26, 0.30), ("west", 0.14, 14, -0.05),
          ("central", 0.18, 23, 0.05)]


def rng(t):
    x = ((t ^ (t >> 16)) * 2246822507) & M32
    x = ((x ^ (x >> 13)) * 3266489909) & M32
    s = [(x ^ (x >> 16)) & M32]

    def next_value():
        s[0] = (s[0] * 1664525 + 1013904223) & M32
        return s[0] / 4294967296
    return next_value


def pick(u, w):
    t = 0.0
    for i, x in enumerate(w):
        t += x
        if u < t:
            return i
    return len(w) - 1


SUBSCRIBERS, DATA = [], []
for i in range(RAW):
    r = rng(SEED + i)
    b = REGION[pick(r(), [x[1] for x in REGION])]
    household = 1 + pick(r(), [0.06, 0.24, 0.30, 0.24, 0.11, 0.05])
    rating = 1 + pick(r(), [0.08, 0.14, 0.27, 0.34, 0.17])
    if r() >= 0.046:
        SUBSCRIBERS.append({"id": 10001 + i, "region": b, "household": household,
                            "satisfaction": rating})
for k in SUBSCRIBERS:
    r, v = rng(SEED + 7000 + k["id"]), []
    for n in range(pick(r(), [0.05, 0.12, 0.21, 0.62])):
        v.append(0.0 if r() < 0.038 else math.floor(
            k["region"][2] * math.exp((r() + r() + r() - 1.5) * 0.62)
            * (1 - n * 0.05) * 100 + 0.5) / 100)
        r()
    if not v:
        continue
    mean, vol, r = sum(v) / len(v), max(v) - min(v), rng(SEED + 51000 + k["id"])
    z = (0.052 * (mean - 20) + 0.85 * (min(v) == 0.0) + 0.026 * vol
         + 0.24 * (k["household"] >= 5) - 0.20 * (k["satisfaction"] >= 4) + k["region"][3]
         + 0.9 * (k["id"] - 10001) / 1399 + (r() + r() + r() - 1.5) * 1.30)
    x = {"household": k["household"], "satisfaction": k["satisfaction"],
         "period": len(v), "mean_consumption": round(mean, 2), "volatility": round(vol, 2),
         "region_name": k["region"][0],
         "zero_readings": int(min(v) == 0.0), "suspicious": int(z > 1.35)}
    for b in REGION:
        x["area_" + b[0]] = int(k["region"][0] == b[0])
    x["per_person"] = round(x["mean_consumption"] / x["household"], 3)
    DATA.append(x)


NUMERIC = ["mean_consumption", "volatility", "household", "satisfaction", "period", "per_person"]
FIELDS = NUMERIC + ["zero_readings", "area_north", "area_south", "area_east", "area_west",
                    "area_central"]


def shuffle(data, seed):
    r, s = rng(seed), list(range(len(data)))
    for i in range(len(s) - 1, 0, -1):
        j = int(r() * (i + 1))
        s[i], s[j] = s[j], s[i]
    return [data[i] for i in s]


K = shuffle(DATA, SEED + 90000)
TRAIN, VAL, TEST = K[:756], K[756:1008], K[1008:]
BASELINE_TEST = sum(x["suspicious"] == 0 for x in TEST) / len(TEST)
SCALE = {a: (statistics.fmean(x[a] for x in TRAIN),
             statistics.pstdev([x[a] for x in TRAIN])) for a in NUMERIC}


def vectorize(s):                           # labels +1 / -1
    return [([(x[a] - SCALE[a][0]) / SCALE[a][1] if a in SCALE else x[a]
              for a in FIELDS], 1 if x["suspicious"] else -1) for x in s]


V_TRAIN, V_VAL, V_TEST = vectorize(TRAIN), vectorize(VAL), vectorize(TEST)


def score(w, b, v):
    return sum(w[i] * v[i] for i in range(len(w))) + b


def train_margin(train, lam, rounds=300, n=None):   # MODEL: margin-based classifier
    d, n = len(train[0][0]), n or len(train)
    w, b = [0.0] * d, 0.0
    for t in range(1, rounds + 1):
        step, gw, gb = 1.0 / (lam * t), [0.0] * d, 0.0
        for v, y in train:
            if y * score(w, b, v) < 1:      # inside the margin or on the wrong side
                for i in range(d):
                    gw[i] -= y * v[i]
                gb -= y
        for i in range(d):
            w[i] -= step * (gw[i] / n + lam * w[i])
        b -= step * gb / n
    return w, b


def accuracy(m, target):
    return sum((score(m[0], m[1], v) > 0) == (y > 0) for v, y in target) / len(target)


CANDIDATES = {}
print(f"baseline on test {BASELINE_TEST:.4f}; training {len(V_TRAIN)} rows, {len(FIELDS)} columns")
print(f"{'penalty':>8} {'|w|':>7} {'margin':>8} {'support':>7} {'validation':>10} {'test':>8}")
for lam in (0.0003, 0.001, 0.003, 0.01, 0.03, 0.1, 0.3, 1.0):
    m = train_margin(V_TRAIN, lam)
    nw = math.sqrt(sum(x * x for x in m[0]))
    CANDIDATES[lam] = (accuracy(m, V_VAL), accuracy(m, V_TEST))
    inside = sum(1 for v, y in V_TRAIN if y * score(m[0], m[1], v) < 1)
    print(f"{lam:>8} {nw:>7.3f} {2 / nw:>8.3f} {inside:>7} "
          f"{CANDIDATES[lam][0]:>10.4f} {CANDIDATES[lam][1]:>8.4f}")
best = max(CANDIDATES, key=lambda a: (CANDIDATES[a][0], a))
smallest = max(CANDIDATES, key=lambda a: (CANDIDATES[a][0], -a))
print(f"candidates tried {len(CANDIDATES)}; penalty chosen on validation {best} -> test "
      f"{CANDIDATES[best][1]:.4f}, {CANDIDATES[best][1] - BASELINE_TEST:+.4f} above the baseline")
print(f"had the tie been resolved in favor of the smallest penalty, penalty {smallest} -> test "
      f"{CANDIDATES[smallest][1]:.4f}")

baseline on test 0.7579; training 756 rows, 12 columns
 penalty     |w|   margin support validation     test
  0.0003  42.768    0.047     252     0.7421   0.7619
   0.001  12.742    0.157     338     0.8294   0.8016
   0.003   4.289    0.466     320     0.8294   0.7976
    0.01   1.587    1.260     324     0.8294   0.7857
    0.03   0.833    2.401     333     0.8175   0.7857
     0.1   0.266    7.523     334     0.8135   0.7579
     0.3   0.075   26.789     335     0.8135   0.7579
     1.0   0.021   93.488     335     0.8135   0.7579
candidates tried 8; penalty chosen on validation 0.01 -> test 0.7857, +0.0278 above the baseline
had the tie been resolved in favor of the smallest penalty, penalty 0.001 -> test 0.8016

The margin widens from 0.047 to 93.488 and both ends of the table are bad. At the narrowest margin the length of the weight vector climbs to 42.768 and the validation number falls 0.0714 below its own baseline; at the widest margin the boundary flattens and on the test set the model reads exactly the baseline. The three candidates in between are tied on the validation set at 0.8294.

Resolving the tie is a decision and the last two lines give its cost. When it is resolved in favor of the wide margin, penalty 0.01 is chosen and the family number is 0.7857, 0.0278 above the baseline; had it been resolved in favor of the narrow margin, 0.8016 would have been reported. The 0.0159 between them is four subscribers and comes from no measurement, only from the tie rule.

How Many Examples Carry the Decision

The number of support vectors is in the table but it is not a claim on its own; whether they carry the decision is measured by dropping those rows and retraining the model.

CHOSEN = 0.01
M = train_margin(V_TRAIN, CHOSEN)
P0 = [score(M[0], M[1], v) > 0 for v, _ in V_TEST]
INSIDE = [i for i, (v, y) in enumerate(V_TRAIN) if y * score(M[0], M[1], v) < 1]
OUTSIDE = [i for i, (v, y) in enumerate(V_TRAIN) if y * score(M[0], M[1], v) >= 1]
r = rng(SEED + 610000)


def random_pick(pool, k):
    h = list(pool)
    for i in range(len(h) - 1, 0, -1):
        j = int(r() * (i + 1))
        h[i], h[j] = h[j], h[i]
    return set(h[:k])


print(f"penalty {CHOSEN}: {len(INSIDE)} rows inside the margin, {len(OUTSIDE)} rows outside the margin")
print(f"{'dropped':>7} {'pool':>14} {'changed predictions':>19} {'test':>8}")
for k in (200, 300):
    for name, pool in (("outside margin", OUTSIDE), ("inside margin", INSIDE)):
        dropped = random_pick(pool, k)
        m2 = train_margin([p for i, p in enumerate(V_TRAIN) if i not in dropped], CHOSEN,
                          n=len(V_TRAIN))
        p2 = [score(m2[0], m2[1], v) > 0 for v, _ in V_TEST]
        print(f"{k:>7} {name:>14} {sum(a != b for a, b in zip(P0, p2)):>19} "
              f"{accuracy(m2, V_TEST):>8.4f}")

penalty 0.01: 324 rows inside the margin, 432 rows outside the margin
dropped           pool changed predictions     test
    200 outside margin                   0   0.7857
    200  inside margin                  15   0.7817
    300 outside margin                   4   0.7937
    300  inside margin                  21   0.7579

In the chosen candidate, 324 of the 756 training rows are inside the margin; this is a measure of how intertwined the classes are. When 300 rows are dropped from outside the margin, the prediction changes on only four of the 252 test rows, and when 200 rows are dropped none changes. When the same number of rows is dropped from inside the margin, the number of changed predictions rises to 21. The number of changed predictions is not a loss measure either: when 300 rows are dropped from outside the margin, the test number rises to 0.7937. What is measured is the movement of the boundary.

The Kernel Transform

A curved boundary is obtained in two ways. The first is opening columns: squares and products are computed and added to the table. The second is the kernel transform: since the model uses nothing but the inner product of two rows in the expanded space, writing a function that gives that inner product directly makes building the space unnecessary. For this the model is written in the dual form: instead of a weight vector, the coefficients of the training rows are kept.

def gram(A, B):                             # inner product matrix; kernels derive from it
    return [[sum(p * q for p, q in zip(a, b)) for b, _ in B] for a, _ in A]


def kernel(name, G, n1, n2):
    if name == "linear":
        return G
    if name == "polynomial":                # (x . z + 1)^2 / 10
        return [[(g + 1.0) ** 2 / 10.0 for g in s] for s in G]
    return [[math.exp(-0.05 * (n1[i] + n2[j] - 2 * g)) for j, g in enumerate(s)]
            for i, s in enumerate(G)]       # radial kernel


def train_dual(K, y, lam, rounds=40):       # MODEL: the same criterion, in dual form
    n, a, b = len(y), [0.0] * len(y), 0.0
    for t in range(1, rounds + 1):
        step = 1.0 / (lam * t)
        f = [sum(a[j] * K[i][j] for j in range(n)) + b for i in range(n)]
        a, gb = [x * (1.0 - step * lam) for x in a], 0.0
        for i in range(n):
            if y[i] * f[i] < 1:
                a[i] += step * y[i] / n
                gb += y[i]
        b += step * gb / n
    return a, b


def dual_accuracy(a, b, KH, target):
    return sum((sum(a[j] * KH[i][j] for j in range(len(a))) + b > 0) == (y > 0)
               for i, (v, y) in enumerate(target)) / len(target)


Y = [y for _, y in V_TRAIN]
G_TRAIN, G_VAL, G_TEST = gram(V_TRAIN, V_TRAIN), gram(V_VAL, V_TRAIN), gram(V_TEST, V_TRAIN)
N_TRAIN = [sum(p * p for p in v) for v, _ in V_TRAIN]
N_VAL = [sum(p * p for p in v) for v, _ in V_VAL]
N_TEST = [sum(p * p for p in v) for v, _ in V_TEST]
KERNELS = {}
print(f"{'kernel':<12}{'penalty':>8}{'validation':>11}{'test':>8}")
for name in ("linear", "polynomial", "radial"):
    K = kernel(name, G_TRAIN, N_TRAIN, N_TRAIN)
    KV, KT = kernel(name, G_VAL, N_VAL, N_TRAIN), kernel(name, G_TEST, N_TEST, N_TRAIN)
    for lam in (0.003, 0.01, 0.03):
        a, b = train_dual(K, Y, lam)
        KERNELS[(name, lam)] = (dual_accuracy(a, b, KV, V_VAL), dual_accuracy(a, b, KT, V_TEST))
        print(f"{name:<12}{lam:>8}{KERNELS[(name, lam)][0]:>11.4f}{KERNELS[(name, lam)][1]:>8.4f}")
chosen = max(KERNELS, key=lambda a: (KERNELS[a][0], -a[1]))
print(f"candidates tried {len(CANDIDATES) + len(KERNELS)} (linear {len(CANDIDATES)} + kernel "
      f"{len(KERNELS)}); chosen on validation {chosen[0]}, penalty {chosen[1]} -> test "
      f"{KERNELS[chosen][1]:.4f}")
print(f"kernel matrix {len(V_TRAIN)} x {len(V_TRAIN)} = {len(V_TRAIN) ** 2} entries")

kernel       penalty validation    test
linear         0.003     0.8016  0.7817
linear          0.01     0.8095  0.7897
linear          0.03     0.8294  0.7857
polynomial     0.003     0.7778  0.7857
polynomial      0.01     0.7897  0.7778
polynomial      0.03     0.8175  0.7937
radial         0.003     0.7222  0.6627
radial          0.01     0.8135  0.7579
radial          0.03     0.8135  0.7579
candidates tried 17 (linear 8 + kernel 9); chosen on validation linear, penalty 0.03 -> test 0.7857
kernel matrix 756 x 756 = 571536 entries

None of the nine candidates can beat the straight boundary. The choice read from validation is the linear kernel and the test number is 0.7857, the same number the primal form gave. The polynomial kernel’s highest test number is 0.7937 and it exceeds the linear kernel’s, but its validation number stays at 0.8175; choosing that candidate would be looking at the test set. The radial kernel at the smallest penalty falls 0.0952 below the baseline. The difference the curved boundary makes on this set is zero when the choice is read from the validation set; its cost is visible, though: the dual form keeps 756 coefficients and computes 571,536 kernel entries.

What remains is the kernel’s real claim: the same model should give the same decision when the columns are computed explicitly.

def explicit_polynomial(v):  # explicit columns equivalent to the (x . z + 1)^2 / 10 kernel
    d = len(v)
    o = [math.sqrt(2) * v[i] for i in range(d)] + [v[i] * v[i] for i in range(d)]
    o += [math.sqrt(2) * v[i] * v[j] for i in range(d) for j in range(i + 1, d)]
    return [x / math.sqrt(10.0) for x in o]


X_TRAIN = [(explicit_polynomial(v), y) for v, y in V_TRAIN]
X_TEST = [(explicit_polynomial(v), y) for v, y in V_TEST]
LAM = 0.03
ma = train_margin(X_TRAIN, LAM, rounds=40)
KT_, KS_ = kernel("polynomial", G_TRAIN, N_TRAIN, N_TRAIN), kernel("polynomial", G_TEST, N_TEST, N_TRAIN)
a, b = train_dual(KT_, Y, LAM)
tk = [sum(a[j] * KS_[i][j] for j in range(len(a))) + b > 0 for i in range(len(V_TEST))]
ta = [score(ma[0], ma[1], v) > 0 for v, _ in X_TEST]
print(f"explicitly computed columns {len(X_TRAIN[0][0])}, columns opened by the kernel 0")
print(f"test with explicit columns  {accuracy(ma, X_TEST):.4f}")
print(f"test with the kernel        {KERNELS[('polynomial', LAM)][1]:.4f}")
print(f"test rows on which the two forms give the same prediction "
      f"{sum(p == q for p, q in zip(tk, ta))}/{len(V_TEST)}")

explicitly computed columns 90, columns opened by the kernel 0
test with explicit columns  0.7937
test with the kernel        0.7937
test rows on which the two forms give the same prediction 252/252

All the second-degree terms of twelve columns open 90 columns and the two forms give the same prediction on 252 of the 252 test rows. The gain is not in the number but in the columns that are not opened: had the degree gone up to three, the number of columns to open would go to 454 while the kernel function would remain a single change of exponent. Seventeen candidates were tried, the family number is 0.7857, and the curved boundary added nothing to the difference above the baseline.

Summary

• The margin width goes from 0.047 to 93.488 with the penalty coefficient; at the narrowest margin the validation number falls 0.0714 below its own baseline.
• The family number is 0.7857 and it beats the baseline by 0.0278; because three candidates tied on validation, 0.0159 of this number comes from the tie rule.
• The decision is carried by the 324 rows inside the margin: when 300 rows are dropped from outside the margin 4 predictions change, from inside 21.
• Nine candidates were tried with three kernels; the chosen kernel is linear and the curved boundary added nothing.
• The polynomial kernel gives the same prediction on 252 of 252 test rows as 90 explicitly computed columns; the cost of the dual form is 571,536 kernel entries.

Next Step

The three families shared one constraint: the sum of distances, the product of probabilities and the inner product all use every column at once, and which column carries the decision can be read only indirectly. The next lesson turns to a family that handles the columns one at a time: the decision tree chooses a single column and a single threshold at each step and makes the choice according to a split criterion. The tree itself was used and its depth measured throughout the Introduction to Machine Learning and Data Preparation and Feature Engineering courses; what has not been asked is which criterion chooses the split and what pruning a grown tree takes back, instead of limiting the depth in advance.