|
| 1 | +#!/usr/bin/env python3 |
| 2 | + |
| 3 | +import collections |
| 4 | +import json |
| 5 | + |
| 6 | +import matplotlib.pyplot as plt |
| 7 | +import numpy as np |
| 8 | +from scipy.optimize import minimize |
| 9 | +from scipy.special import logit |
| 10 | + |
| 11 | +from codeclash import REPO_DIR |
| 12 | + |
| 13 | +data = json.loads((REPO_DIR / "round_scores.json").read_text()) |
| 14 | +onlyrows = collections.defaultdict(list) |
| 15 | +for game in data.keys(): |
| 16 | + for model_results in data[game]: |
| 17 | + for model in model_results.keys(): |
| 18 | + if len(model_results[model]) != 15: |
| 19 | + # print("skipping", game, len(model_results[model])) |
| 20 | + continue |
| 21 | + onlyrows[game].append(model_results[model]) |
| 22 | + |
| 23 | + |
| 24 | +onlyrows = {game: np.array(onlyrows[game]) for game in onlyrows.keys()} |
| 25 | +onlyrows["BattleSnake"].shape |
| 26 | +plt.hist(onlyrows["BattleSnake"][:, 14]) |
| 27 | + |
| 28 | + |
| 29 | +METHOD = "ar1" |
| 30 | +# METHOD = "arp_geometric_decay" |
| 31 | +# METHOD = "ar_p" |
| 32 | +# METHOD = "same_as_previous" |
| 33 | + |
| 34 | + |
| 35 | +def ar_residuals_logit(scores: np.ndarray, gammas: np.ndarray, mus: np.ndarray) -> np.ndarray: |
| 36 | + """ |
| 37 | + One-step-ahead residuals for AR(p) on logit(scores) with stationary means. |
| 38 | +
|
| 39 | + scores: (T, R) in (0,1) |
| 40 | + gammas: (p,) |
| 41 | + mus: (T,) |
| 42 | + Returns: residuals (T, R) with NaN for the first p columns. |
| 43 | + """ |
| 44 | + eps = 1e-6 |
| 45 | + Z = logit(np.clip(scores, eps, 1 - eps)) |
| 46 | + gammas = np.atleast_1d(gammas).astype(float) |
| 47 | + p = len(gammas) |
| 48 | + T, R = Z.shape |
| 49 | + res = np.full((T, R), np.nan, dtype=float) |
| 50 | + mu = mus[:, None] |
| 51 | + for t in range(p, R): |
| 52 | + hat = mu.copy() |
| 53 | + for k in range(1, p + 1): |
| 54 | + hat += gammas[k - 1] * (Z[:, t - k : t - k + 1] - mu) |
| 55 | + res[:, t] = (Z[:, t : t + 1] - hat).ravel() |
| 56 | + return res |
| 57 | + |
| 58 | + |
| 59 | +def fit(scores: np.ndarray, p: int = 1): |
| 60 | + """ |
| 61 | + Fit AR(p) on logit(scores) by conditional Gaussian MLE (CSS). |
| 62 | +
|
| 63 | + scores: (T, R) |
| 64 | + p: AR order (default 1) |
| 65 | + Returns: (gammas_hat: (p,), mus_hat: (T,), residual_std: float, result) |
| 66 | + """ |
| 67 | + T, R = scores.shape |
| 68 | + if R <= p: |
| 69 | + raise ValueError("Number of rounds must exceed p.") |
| 70 | + |
| 71 | + eps = 1e-6 |
| 72 | + Z = logit(np.clip(scores, eps, 1 - eps)) |
| 73 | + |
| 74 | + mu0 = Z.mean(axis=1) # (T,) |
| 75 | + game0 = np.full(p, 0.3, dtype=float) |
| 76 | + |
| 77 | + def pack(game: np.ndarray, mu: np.ndarray) -> np.ndarray: |
| 78 | + return np.concatenate([game, mu]) |
| 79 | + |
| 80 | + def unpack(theta: np.ndarray) -> tuple[np.ndarray, np.ndarray]: |
| 81 | + game = theta[:p] |
| 82 | + mu = theta[p:] |
| 83 | + return game, mu |
| 84 | + |
| 85 | + valid_mask = np.ones((T, R), dtype=bool) |
| 86 | + valid_mask[:, :p] = False |
| 87 | + N = valid_mask.sum() |
| 88 | + |
| 89 | + def nll(theta: np.ndarray) -> float: |
| 90 | + game, mu = unpack(theta) |
| 91 | + res = ar_residuals_logit(scores, game, mu) |
| 92 | + e = res[valid_mask] |
| 93 | + sigma2 = (e @ e) / N |
| 94 | + return 0.5 * N * np.log(sigma2) |
| 95 | + |
| 96 | + theta0 = pack(game0, mu0) |
| 97 | + lower_game = np.full(p, -0.999) |
| 98 | + upper_game = np.full(p, 0.999) |
| 99 | + mu_bounds = [(-np.inf, np.inf)] * T |
| 100 | + bounds = list(zip(lower_game, upper_game)) + mu_bounds |
| 101 | + |
| 102 | + result = minimize(nll, theta0, method="L-BFGS-B", bounds=bounds) |
| 103 | + game_hat, mu_hat = unpack(result.x) |
| 104 | + res_hat = ar_residuals_logit(scores, game_hat, mu_hat) |
| 105 | + residual_std = np.nanstd(res_hat) |
| 106 | + |
| 107 | + return game_hat, mu_hat, residual_std, result |
| 108 | + |
| 109 | + |
| 110 | +def residuals_same_as_previous(scores: np.ndarray): |
| 111 | + """ |
| 112 | + Compute residuals for a model where each tournament is predicted using |
| 113 | + the logit of the previous round's score. |
| 114 | +
|
| 115 | + scores: shape (num_tournaments, num_rounds), entries in (0,1) |
| 116 | + Returns: residuals array of shape (num_tournaments, num_rounds) |
| 117 | + """ |
| 118 | + eps = 1e-6 |
| 119 | + Z = logit(np.clip(scores, eps, 1 - eps)) |
| 120 | + |
| 121 | + num_tournaments, num_rounds = Z.shape |
| 122 | + residuals = np.zeros_like(Z) |
| 123 | + residuals[:, 0] = 0 |
| 124 | + residuals[:, 1:] = Z[:, 1:] - Z[:, :-1] |
| 125 | + return residuals |
| 126 | + |
| 127 | + |
| 128 | +def residuals_same_as_first(scores: np.ndarray): |
| 129 | + """ |
| 130 | + Compute residuals for a model where each tournament is predicted using |
| 131 | + the logit of the first round's score. |
| 132 | +
|
| 133 | + scores: shape (num_tournaments, num_rounds), entries in (0,1) |
| 134 | + Returns: residuals array of shape (num_tournaments, num_rounds) |
| 135 | + """ |
| 136 | + eps = 1e-6 |
| 137 | + Z = logit(np.clip(scores, eps, 1 - eps)) |
| 138 | + return Z - Z[:, [0]] |
| 139 | + |
| 140 | + |
| 141 | +def residuals_all_5050(scores: np.ndarray): |
| 142 | + """ |
| 143 | + Compute residuals for a model where each tournament is predicted to be 50% 50% |
| 144 | + """ |
| 145 | + num_tournaments, num_rounds = scores.shape |
| 146 | + eps = 1e-6 |
| 147 | + Z = logit(np.clip(scores, eps, 1 - eps)) |
| 148 | + return Z - logit(0.5) |
| 149 | + |
| 150 | + |
| 151 | +def residuals_all_mean(scores: np.ndarray): |
| 152 | + """ |
| 153 | + Compute residuals for a model where each tournament is predicted to be the mean of the scores |
| 154 | + """ |
| 155 | + eps = 1e-6 |
| 156 | + Z = logit(np.clip(scores, eps, 1 - eps)) |
| 157 | + mean_logits = Z.mean(axis=1, keepdims=True) |
| 158 | + return Z - mean_logits |
| 159 | + |
| 160 | + |
| 161 | +print("=== baselines ===") |
| 162 | +print("same_as_previous") |
| 163 | +print(residuals_same_as_previous(onlyrows["BattleSnake"]).std()) |
| 164 | +print("same_as_first") |
| 165 | +print(residuals_same_as_first(onlyrows["BattleSnake"]).std()) |
| 166 | +print("all_5050") |
| 167 | +print(residuals_all_5050(onlyrows["BattleSnake"]).std()) |
| 168 | +print("all_mean") |
| 169 | +print(residuals_all_mean(onlyrows["BattleSnake"]).std()) |
| 170 | + |
| 171 | +print("=== fitted ===") |
| 172 | + |
| 173 | +print("p=1") |
| 174 | +game_hat, mu_hat, residual_std, result = fit(onlyrows["BattleSnake"], p=1) |
| 175 | +print(game_hat, residual_std) |
| 176 | + |
| 177 | +print("p=2") |
| 178 | +game_hat, mu_hat, residual_std, result = fit(onlyrows["BattleSnake"], p=2) |
| 179 | +print(game_hat, residual_std) |
| 180 | + |
| 181 | +print("p=3") |
| 182 | +game_hat, mu_hat, residual_std, result = fit(onlyrows["BattleSnake"], p=3) |
| 183 | +print(game_hat, residual_std) |
0 commit comments