Hyperparameter Search

The previous lesson measured a single model, and the number of candidates tried was one. The rule of the course is this: the number of a model cannot be read without writing how many candidates were tried to find that number; a number whose tuning budget is not written counts as unmeasured. This lesson is where that rule is paid.

The measure has two parts. The first is the search itself: when the same space and the same budget are given, which candidate do the three procedures — grid, random and model-based — find. The second is what the budget inflates: as the number of candidates rises from 1 to 200, what opens between the best number on validation and the test number of that candidate. This inflation is finally put side by side with the difference added by changing the family.

• ET50. The tuning space has four axes: the penalty coefficient ($10^{-3}$–$10^{3}$, continuous on a log scale), the decision threshold (0.20–0.60, continuous), the polynomial degree (1 or 2) and the number of columns taken from a fixed ordering (3, 6 or 9). All three searches work in the same space.
• ET51. The choice of candidate is read only from the validation set (252 rows); the test set enters no choice and is read only to report the number of the selected candidate.
• ET52. The budget unit is the number of candidates trained and measured on the validation set; the proposals that model-based search tries on the surrogate model need no training and are not counted. Because the product of the discrete axes (2 degrees × 3 columns = 6) divides the budget, grid search cannot fill the given budget exactly.
• ET53. Model-based search builds a surrogate model from the validation numbers of the candidates already tried: a distance-weighted mean plus an exploration bonus proportional to the distance to the nearest candidate tried. This is a model too.
• ET54. The second family is the single-split decision rule: one column, one threshold, one direction. It is structurally weaker than the linear model and is searched with the same budget.
• ET55. The baseline carries over from the Introduction to Machine Learning course: the procedure that says the most frequent class, 0.7579 on the test set. The number of candidates is printed at the end of the lesson and repeated candidates are counted once.

Same Space, Same Budget, Three Searches

# COMMON — MODEL: The FICTIONAL subscriber table of the Introduction to Machine Learning course
# is rebuilt with the same seed; the model is a linear classifier fitted to the 0/1 label by
# least squares.
import math

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 d 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 - d * 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 = {"id": k["id"], "household": k["household"], "satisfaction": k["satisfaction"],
         "period": len(v), "mean_consumption": round(mean, 2), "volatility": round(vol, 2),
         "zero_readings": int(min(v) == 0.0), "suspicious": int(z > 1.35),
         "region_name": k["region"][0]}
    for b in REGION:
        x["area_" + b[0]] = int(k["region"][0] == b[0])
    DATA.append(x)


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:]
ORDER = ["mean_consumption", "zero_readings", "volatility", "household", "satisfaction",
         "period", "area_east", "area_west", "area_south"]
SQUARED = {"mean_consumption", "volatility", "household", "satisfaction", "period"}
COLUMNS = (3, 6, 9)


def select_fields(j, degree):
    a = ORDER[:j]
    return a + ([q + "^2" for q in a if q in SQUARED] if degree == 2 else [])


def value(x, name):
    return float(x[name[:-2]]) ** 2 if name.endswith("^2") else float(x[name])


def solve(A, b):
    n = len(A)
    M = [list(A[i]) + [b[i]] for i in range(n)]
    for i in range(n):
        p = max(range(i, n), key=lambda r: abs(M[r][i]))
        M[i], M[p] = M[p], M[i]
        for r in range(n):
            if r != i:
                f = M[r][i] / M[i][i]
                for c in range(i, n + 1):
                    M[r][c] -= f * M[i][c]
    return [M[i][n] / M[i][i] for i in range(n)]


def fit(s, fields, penalty):
    d = len(fields) + 1
    A = [[0.0] * d for _ in range(d)]
    b = [0.0] * d
    for x in s:
        v = [1.0] + [value(x, a) for a in fields]
        for i in range(d):
            b[i] += v[i] * x["suspicious"]
            for j in range(i, d):
                A[i][j] += v[i] * v[j]
    for i in range(d):
        for j in range(i):
            A[i][j] = A[j][i]
        if i:
            A[i][i] += penalty
    return solve(A, b)


def accuracy(w, fields, threshold, s):
    d = 0
    for x in s:
        p = w[0] + sum(w[i + 1] * value(x, a) for i, a in enumerate(fields))
        d += (1 if p >= threshold else 0) == x["suspicious"]
    return d / len(s)


CACHE, TRIED = {}, set()


def evaluate(candidate):
    penalty, threshold, degree, j = candidate
    k = (round(penalty, 9), degree, j)
    if k not in CACHE:
        CACHE[k] = (fit(TRAIN, select_fields(j, degree), penalty), select_fields(j, degree))
    TRIED.add(("linear", k, round(threshold, 9)))
    w, fields = CACHE[k]
    return accuracy(w, fields, threshold, VAL), accuracy(w, fields, threshold, TEST)


def best_of(candidates, measure):
    best = None
    for a in candidates:
        d, s = measure(a)
        if best is None or d > best[0]:
            best = (d, s, a)
    return best


c = int(sum(x["suspicious"] for x in TRAIN) * 2 > len(TRAIN))
BASELINE = sum(x["suspicious"] == c for x in TEST) / len(TEST)
DEFAULT_VAL, DEFAULT_TEST = evaluate((1.0, 0.5, 1, 9))
print(f"baseline {BASELINE:.4f}")
print(f"untuned default candidate: validation {DEFAULT_VAL:.4f}, test {DEFAULT_TEST:.4f}")

baseline 0.7579
untuned default candidate: validation 0.8452, test 0.7976

The untuned default candidate — penalty 1, threshold 0.5, degree 1, nine columns — gives 0.8452 on the validation set. The three searches below and the 200-candidate sweep will try to beat this number.

def random_candidate(r):
    return (10 ** (-3 + 6 * r()), 0.20 + 0.40 * r(),
            1 + int(r() * 2), COLUMNS[int(r() * 3)])


def grid(budget):
    top = (1, 1)
    for n_pen in range(1, budget + 1):
        for n_thr in range(1, budget + 1):
            if n_pen * n_thr * 6 <= budget and \
                    (n_pen * n_thr, -abs(n_pen - n_thr), n_pen) > \
                    (top[0] * top[1], -abs(top[0] - top[1]), top[0]):
                top = (n_pen, n_thr)
    n_pen, n_thr = top
    return [(10 ** (-3 + 6 * (i + 0.5) / n_pen), 0.20 + 0.40 * (t + 0.5) / n_thr, degree, j)
            for i in range(n_pen) for t in range(n_thr) for degree in (1, 2) for j in COLUMNS]


def random_search(budget, seed=SEED + 4100):
    r = rng(seed)
    return [random_candidate(r) for _ in range(budget)]


def distance(a, b):
    return (((math.log10(a[0]) - math.log10(b[0])) / 6) ** 2
            + ((a[1] - b[1]) / 0.4) ** 2 + (a[2] - b[2]) ** 2
            + ((a[3] - b[3]) / 6) ** 2)


def model_based(budget, seed=SEED + 4200, initial=8, proposals=40):
    r, history = rng(seed), []
    for _ in range(min(initial, budget)):
        a = random_candidate(r)
        history.append((a, evaluate(a)[0]))
    while len(history) < budget:
        top_score, top_candidate = None, None
        for _ in range(proposals):
            a = random_candidate(r)
            w = [1.0 / (distance(a, g) + 0.01) for g, _ in history]
            score = (sum(q * s for q, (_, s) in zip(w, history)) / sum(w)
                     + 0.02 * min(distance(a, g) for g, _ in history) ** 0.5)
            if top_score is None or score > top_score:
                top_score, top_candidate = score, a
        history.append((top_candidate, evaluate(top_candidate)[0]))
    return [g for g, _ in history]


print(f"{'search':<16} {'candidates':>10} {'validation':>10} {'test':>8} {'inflation':>9}  setting")
for name, fn in (("grid", grid), ("random", random_search),
                 ("model-based", model_based)):
    A = fn(50)
    d, s, a = best_of(A, evaluate)
    print(f"{name:<16} {len(A):>10} {d:>10.4f} {s:>8.4f} {d - s:>+9.4f}  "
          f"penalty {a[0]:.3g}, threshold {a[1]:.3f}, degree {a[2]}, columns {a[3]}")

print(f"\n{'budget':>6} {'validation':>10} {'test':>8} {'inflation':>9}")
R = random_search(200)
for b in (1, 2, 5, 10, 20, 50, 100, 200):
    d, s, a = best_of(R[:b], evaluate)
    print(f"{b:>6} {d:>10.4f} {s:>8.4f} {d - s:>+9.4f}")

search           candidates validation     test inflation  setting
grid                     48     0.8452   0.7937   +0.0516  penalty 5.62, threshold 0.500, degree 2, columns 9
random                   50     0.8373   0.7857   +0.0516  penalty 0.00481, threshold 0.495, degree 1, columns 9
model-based              50     0.8452   0.8016   +0.0437  penalty 0.00103, threshold 0.540, degree 1, columns 9

budget validation     test inflation
     1     0.7937   0.7778   +0.0159
     2     0.7937   0.7778   +0.0159
     5     0.8254   0.8056   +0.0198
    10     0.8294   0.7976   +0.0317
    20     0.8294   0.7976   +0.0317
    50     0.8373   0.7857   +0.0516
   100     0.8373   0.7857   +0.0516
   200     0.8452   0.7976   +0.0476

In the first table grid search can use 48 candidates against a budget of 50. The product of the discrete axes is six, only 8 points are left for the continuous axes, and these split 4×2: four values on the penalty axis, that is, six orders of magnitude are represented by four points. Random search with the same budget tries 50 distinct values on that axis.

The test numbers of the three searches are 0.7937, 0.7857 and 0.8016 — not the same as the validation ranking. Model-based search is ahead on both sets and has the smallest inflation (+0.0437), but the test difference among the three is 0.0159, that is, four rows; there is not enough evidence to attribute a difference of this size to the nature of the searches.

The main finding is in the second table. The validation column rises by definition: as the budget grows, the best number is the maximum of a growing set and cannot shrink. From 0.7937 to 0.8452, a total of +0.0515. The test column does not follow this climb: from 0.7778 to 0.7976, a total of +0.0198. 62% of the climb never passes to the test set, and in the inflation column it opens from +0.0159 to +0.0476.

The cause is that the validation set is 252 rows. Each row is worth 0.0040; when the real difference between two candidates is below this resolution, which rows happened to fall into that set decides the winner. When 200 candidates are tried, the candidate with the highest number is not the best one but the one that best fit the noise of the validation set — the form, in the tuning space, of the selection bias measured in the Statistical Analysis and Exploratory Data Analysis course.

The point of comparison is in the opening block: the untuned default candidate already gave 0.8452 on the validation set. 200 candidates could not beat this number, and where they reached it they stopped at the same point (0.7976) on the test set as well. The search produced zero gain in return for 200 trainings; we can say this only once the test set is read.

What the Family Adds and What the Budget Inflates

The second claim of the course turns into a number here: the difference added by changing the family is often smaller than the difference inflated by the tuning budget. The second family is structurally weaker — a single-split decision rule: one column, one threshold, one direction.

THRESHOLD_TABLE = {}
for a in ORDER:
    d = sorted({x[a] for x in TRAIN})
    THRESHOLD_TABLE[a] = (d[1:] if len(d) < 12
                          else [d[int(i * len(d) / 12)] for i in range(1, 12)])


def stump_evaluate(candidate):
    a, v, direction = candidate
    TRIED.add(("stump",) + candidate)

    def predict(x):
        return int((x[a] >= v) == (direction == 1))
    return (sum(predict(x) == x["suspicious"] for x in VAL) / len(VAL),
            sum(predict(x) == x["suspicious"] for x in TEST) / len(TEST))


def stump_random(budget, seed=SEED + 4300):
    r, out = rng(seed), []
    for _ in range(budget):
        a = ORDER[int(r() * len(ORDER))]
        out.append((a, THRESHOLD_TABLE[a][int(r() * len(THRESHOLD_TABLE[a]))], 1 + int(r() * 2)))
    return out


print(f"{'budget':>6} {'stump test':>10} {'linear test':>11} {'family difference':>17} "
      f"{'budget inflation':>16}")
for b in (1, 10, 50, 200):
    gd, gs, ga = best_of(stump_random(b), stump_evaluate)
    dd, ds, da = best_of(R[:b], evaluate)
    print(f"{b:>6} {gs:>10.4f} {ds:>11.4f} {ds - gs:>+17.4f} {dd - ds:>+16.4f}")
gd, gs, ga = best_of(stump_random(200), stump_evaluate)
print(f"\nstump budget 200: single split '{ga[0]}', validation {gd:.4f}, test {gs:.4f}")
print(f"number of candidates tried in this lesson {len(TRIED)}")

budget stump test linear test family difference budget inflation
     1     0.7659      0.7778           +0.0119          +0.0159
    10     0.7381      0.7976           +0.0595          +0.0317
    50     0.8016      0.7857           -0.0159          +0.0516
   200     0.8016      0.7976           -0.0040          +0.0476

stump budget 200: single split 'mean_consumption', validation 0.8254, test 0.8016
number of candidates tried in this lesson 352

The family difference column changes sign across the budget: +0.0119 at budget 1, +0.0595 at budget 10, it turns negative at budget 50 (−0.0159) and stays negative at budget 200 (−0.0040). The structurally stronger family could not beat the structurally weaker family on the test set at two budgets. The 0.0595 at budget 10 comes from the stump happening to get stuck on a bad candidate; there the stump is at 0.7381, below the baseline.

The last column carries the difference inflated by the budget on the same rows: +0.0159, +0.0317, +0.0516, +0.0476. At budget 200, while the family difference is −0.0040, the budget inflation is +0.0476 — the inflation is nearly twelve times the family difference and in the opposite direction. Changing the family of the model did not move the number on this set; what moved it was how many candidates were tried, and that movement was not real.

The rule the stump found rests on a single column: a single split on mean consumption, 0.8016 on the test set. The nine-column, penalized linear models of this lesson could not beat a single threshold.

The number of candidates tried is 352 and it is smaller than the sum of the given budgets. There are 299 distinct candidates in the linear family: one default candidate, 148 candidates of the three searches and 150 candidates added by the budget sweep. In the stump family, on the other hand, there are 53 distinct candidates despite a budget of 200, because the whole tuning space of that family has only 74 points (37 thresholds, two directions). This is why budgets 50 and 200 give the same number: the space is saturated and enlarging the budget does not produce new candidates. The tuning budget of a family cannot be read without writing the size of its space.

The reported number is the test-set number of the candidate selected at budget 200: 0.7976, 0.0397 above the baseline.

Summary

• With the same budget of 50 candidates, grid search can use only 48 (the product of the discrete axes divides the budget); the test numbers of the three searches are 0.7937, 0.7857 and 0.8016.
• As the budget rises from 1 to 200, the best number on the validation set climbs by definition from 0.7937 to 0.8452, the test number of the same candidate goes from 0.7778 to 0.7976, and the inflation opens from +0.0159 to +0.0476.
• The untuned default candidate already gives 0.8452 on the validation set; 200 candidates could not beat this number and the search produced zero gain in return for the budget spent.
• At budgets of 50 and 200 the single-split decision rule beats the linear model on the test set: the family difference is −0.0159 and −0.0040 while the budget inflation is +0.0516 and +0.0476, that is, the inflation is nearly twelve times the family difference.
• 352 distinct candidates were tried in this lesson (299 in the linear family, 53 in the stump family); the budget of the stump family was 200 but the search saturated because its space has only 74 points, and enlarging the budget did not produce new candidates. The reported number is the test number of the candidate selected at budget 200 (0.7976).

Next Step

All the numbers of this lesson are a single ratio: in how many of 252 rows the right decision was made. Which regions, which tariff steps, which household sizes and how many readings the rows in the denominator of that ratio came from was never asked. The next lesson asks the question by splitting the error into these slices, and measures that the number of bad slices found also grows by chance as the number of slices increases.