SLD - Parallel Scaling Law - SLDAgent + o4-mini

Best Run 2 R² = 0.999958

▼

Python

import numpy as np
from scipy.optimize import least_squares

def scaling_law_func(data_points, params):
    """
    Four-parameter scaling law with diminishing-parallel returns:
      loss = b + a * (N/1e9)^(-alpha) * (1 + log2(P))^(-beta)
    where:
      N      = model size (num_params)
      P      = parallel_size
      P_eff  = 1 + log2(P) for diminishing returns
    params = [a, alpha, beta, b]
    """
    X = np.atleast_2d(np.asarray(data_points, dtype=float))
    N = X[:, 0] / 1e9
    P = X[:, 1]
    # effective parallel factor
    eps = 1e-12
    P_eff = 1.0 + np.log2(np.clip(P, eps, None))
    a, alpha, beta, b = params
    # enforce minimal positivity for stability
    a     = max(a, eps)
    alpha = max(alpha, eps)
    beta  = max(beta, eps)
    b     = max(b, 0.0)
    return b + a * (N ** (-alpha)) * (P_eff ** (-beta))


def fit_scaling_law(data_points, loss_values):
    """
    Fit the 4-parameter law
      loss = b + a*(N/1e9)^(-alpha)*(1+log2(P))^(-beta)
    via:
      1) init b near the lowest observed loss
      2) log-linear regression for [a, alpha, beta]
      3) robust non-linear least squares (Huber) with analytic Jacobian
    """
    X = np.atleast_2d(np.asarray(data_points, dtype=float))
    y = np.asarray(loss_values, dtype=float).ravel()
    if X.shape[0] != y.size:
        raise ValueError("data_points and loss_values must match lengths")

    # preprocess features
    N = X[:, 0] / 1e9
    P = X[:, 1]
    eps = 1e-12
    P_eff = 1.0 + np.log2(np.clip(P, eps, None))
    y_min = float(np.min(y))

    # 1) initialize intercept b near the lower envelope
    b0 = max(0.0, 0.9 * y_min)

    # 2) log-linear initialization for a, alpha, beta
    y_shift = y - b0
    # clamp to positive
    y_shift = np.clip(y_shift, eps, None)
    logy   = np.log(y_shift)
    logN   = np.log(np.clip(N,   eps, None))
    logPe  = np.log(P_eff)

    # design matrix: logy ≈ C0 − alpha*logN − beta*logPe
    A = np.column_stack([np.ones_like(logy), -logN, -logPe])
    try:
        C0, alpha0, beta0 = np.linalg.lstsq(A, logy, rcond=None)[0]
        a0     = max(np.exp(C0),     eps)
        alpha0 = max(alpha0,         eps)
        beta0  = max(beta0,          eps)
    except Exception:
        span   = max(np.max(y) - y_min, eps)
        a0, alpha0, beta0 = span, 0.5, 0.2

    init_params = np.array([a0, alpha0, beta0, b0], dtype=float)

    # residuals for least_squares (raw domain)
    def residuals(p):
        return scaling_law_func(X, p) - y

    # analytic Jacobian ∂r/∂p
    def jac(p):
        a, alpha, beta, _ = p
        M = (N ** (-alpha)) * (P_eff ** (-beta))
        da     = M
        dalpha = -a * M * np.log(np.clip(N, eps, None))
        dbeta  = -a * M * np.log(P_eff)
        db     = np.ones_like(M)
        return np.vstack([da, dalpha, dbeta, db]).T

    # bounds: a,alpha,beta ≥ eps; b in [0, y_min]
    lower = [eps, eps, eps, 0.0]
    upper = [np.inf, np.inf, np.inf, y_min]

    # 3) robust fitting with Huber loss
    result = least_squares(
        residuals,
        init_params,
        jac=jac,
        bounds=(lower, upper),
        loss='huber',
        f_scale=1e-3,
        xtol=1e-12,
        ftol=1e-12,
        gtol=1e-12
    )

    p_opt = result.x if result.success else init_params
    # ensure intercept never exceeds the observed minimum
    p_opt[3] = min(p_opt[3], y_min)
    return p_opt

#2 Run 1 R² = 0.999954

▼

Python

import numpy as np
from scipy.optimize import least_squares

# EVOLVE-BLOCK-START
def scaling_law_func(data_points, params):
    """
    Predict loss as:
      loss = b + a * num_params^(−alpha) * parallel_size^(−beta)
    using a stable exp(·) form.
    params = [a, alpha, beta, b] (4 parameters)
    """
    D = np.atleast_2d(np.asarray(data_points, dtype=float))
    n, s = D[:, 0], D[:, 1]
    a, alpha, beta, b = params
    # add small eps inside logs for numerical stability
    return b + a * np.exp(-alpha * np.log(n + 1e-12)
                          - beta  * np.log(s + 1e-12))

def fit_scaling_law(data_points, loss_values):
    """
    Fit the 4-parameter scaling law via:
      1) three baseline offsets b₀ ∈ {0.0, 0.5·min(y), 0.9·min(y)}
      2) log-linear least squares to initialize [loga, α, β]
      3) refine all [a, α, β, b] by Levenberg-Marquardt (least_squares)
         with analytic Jacobian on residuals f(data; p) − y
    Returns optimized params [a, alpha, beta, b].
    """
    D = np.atleast_2d(np.asarray(data_points, dtype=float))
    y = np.ravel(np.asarray(loss_values, dtype=float))
    n, s = D[:, 0], D[:, 1]
    # precompute logs
    ln_n = np.log(n + 1e-12)
    ln_s = np.log(s + 1e-12)
    y_min = y.min()
    eps = 1e-12

    best_cost = np.inf
    best_params = None

    # Multi-start on baseline guesses for b
    for frac in (0.0, 0.5, 0.9):
        b0 = max(0.0, frac * y_min)
        # prepare for log-linear init
        y_adj = np.clip(y - b0, eps, None)
        L     = np.log(y_adj)
        A     = np.vstack([np.ones_like(L), -ln_n, -ln_s]).T
        sol, *_ = np.linalg.lstsq(A, L, rcond=None)
        loga0, alpha0, beta0 = sol
        a0 = np.exp(loga0)
        p0 = np.array([max(a0, eps),
                       max(alpha0, 0.0),
                       max(beta0,  0.0),
                       b0], dtype=float)

        # bounds: a>0, α≥0, β≥0, b≥0
        lb = [eps, 0.0, 0.0, 0.0]
        ub = [np.inf, 10.0, 10.0, y_min]

        # residuals and Jacobian for least_squares
        def residuals(p):
            a, alpha, beta, b = p
            # f = b + a * exp(−α ln n − β ln s)
            exp_term = np.exp(-alpha * ln_n - beta * ln_s)
            return b + a * exp_term - y

        def jac(p):
            a, alpha, beta, b = p
            exp_term = np.exp(-alpha * ln_n - beta * ln_s)
            # df/da, df/dα, df/dβ, df/db
            df_da     = exp_term
            df_dalpha = -a * exp_term * ln_n
            df_dbeta  = -a * exp_term * ln_s
            df_db     = np.ones_like(y)
            # stack columns
            return np.vstack([df_da, df_dalpha, df_dbeta, df_db]).T

        try:
            res = least_squares(residuals,
                                p0,
                                jac=jac,
                                bounds=(lb, ub),
                                xtol=1e-12,
                                ftol=1e-12)
            if res.success and res.cost < best_cost:
                best_cost   = res.cost
                best_params = res.x
        except Exception:
            continue

    # fallback if all starts fail
    if best_params is None:
        return np.array([1e-3, 0.5, 0.1, eps], dtype=float)
    return best_params
# EVOLVE-BLOCK-END

#3 Run 3 R² = 0.999954

▼

Python

import numpy as np
from scipy.optimize import minimize

def scaling_law_func(data_points, params):
    """
    Predict loss using a power‐law in model size and parallelism:
      loss = b + a * (num_params/1e9)^(-alpha) * (parallel_size)^(-beta)
    params: [a, alpha, beta, b]
    """
    X = np.atleast_2d(np.asarray(data_points))
    N = X[:, 0] / 1e9
    P = X[:, 1]
    a, alpha, beta, b = params
    return b + a * (N ** (-alpha)) * (P ** (-beta))


def fit_scaling_law(data_points, loss_values):
    """
    Fit a 4-parameter scaling law:
      loss = b + a * (num_params/1e9)^(-alpha) * (parallel_size)^(-beta)
    Steps:
      1) Estimate offset b0 from lower end of the losses.
      2) Solve a log-linear regression with a small ridge to get ln(a), alpha, beta.
      3) Refine offset b0 from residuals.
      4) Jointly optimize [a, alpha, beta, b] by minimizing MSE with L-BFGS-B.
    """
    X = np.atleast_2d(np.asarray(data_points))
    y = np.ravel(loss_values)
    N = X[:, 0] / 1e9
    P = X[:, 1]

    # 1) initial offset estimate
    y_min, y_max = y.min(), y.max()
    b0 = max(0.0, y_min - 0.1 * (y_max - y_min))

    # 2) log-linear regression for ln(a), alpha, beta
    r = np.clip(y - b0, 1e-8, None)
    lnN, lnP = np.log(N), np.log(P)
    # Design matrix: [1, -ln(N), -ln(P)]
    Z = np.vstack([np.ones_like(r), -lnN, -lnP]).T
    t = np.log(r)
    ridge = 1e-6 * np.eye(3)
    theta = np.linalg.solve(Z.T @ Z + ridge, Z.T @ t)
    ln_a0, alpha0, beta0 = theta
    a0 = np.exp(ln_a0)

    # 3) refine b0 from residuals
    pred0 = a0 * (N ** (-alpha0)) * (P ** (-beta0))
    b0 = max(0.0, float(np.mean(y - pred0)))

    init = np.array([a0, alpha0, beta0, b0], dtype=float)

    # 4) final joint optimization on MSE
    def mse_loss(p):
        pred = scaling_law_func(X, p)
        return np.mean((pred - y) ** 2)

    bounds = [(1e-12, None), (1e-8, None), (1e-8, None), (0.0, None)]
    res = minimize(mse_loss, init, method='L-BFGS-B', bounds=bounds)

    return res.x if res.success else init

#4 Run 4 R² = 0.999954

▼

#5 Run 5 R² = 0.999954

▼

Python

import numpy as np
from scipy.optimize import minimize

def scaling_law_func(data_points, params):
    """
    4‐parameter parallel scaling law:
      loss = b + a * (num_params/1e9)^(-alpha) * parallel_size^(-beta)
    params = [a, alpha, beta, b]
    """
    X = np.atleast_2d(np.asarray(data_points, dtype=float))
    # normalize model size (in billions) for numerical stability
    N = X[:, 0] / 1e9
    s = X[:, 1]
    a, alpha, beta, b = params
    return b + a * (N ** (-alpha)) * (s ** (-beta))


def fit_scaling_law(data_points, loss_values):
    """
    Fit the 4-parameter form
      loss = b + a * (N/1e9)^(-alpha) * s^(-beta)
    by a two‐stage procedure:
      1) grid‐search b in [0, 0.8*min(y)] and solve
         log(y-b) ≈ log(a) - alpha*log(N) - beta*log(s)
         via linear least squares to get initial (a,alpha,beta,b).
      2) refine all four parameters by L-BFGS-B minimizing MSE.
    """
    # prepare data
    X = np.atleast_2d(np.asarray(data_points, dtype=float))
    y = np.asarray(loss_values, dtype=float).ravel()
    # normalized feature
    N = X[:, 0] / 1e9
    s = X[:, 1]
    # logs for linear solve
    logN = np.log(N + 1e-12)
    logs = np.log(s + 1e-12)

    y_min = y.min()
    best_err = np.inf
    best_init = None

    # coarse grid for baseline b
    for b in np.linspace(0.0, 0.8 * y_min, 30):
        y_adj = y - b
        # skip invalid offsets
        if np.any(y_adj <= 0):
            continue
        L = np.log(y_adj)
        # design matrix: L = c0 - alpha*logN - beta*logs
        A = np.vstack([np.ones_like(L), -logN, -logs]).T
        # solve for [c0, alpha, beta]
        coeffs, *_ = np.linalg.lstsq(A, L, rcond=None)
        c0, alpha, beta = coeffs
        a = np.exp(c0)
        # compute MSE in original space
        pred = b + a * (N ** (-alpha)) * (s ** (-beta))
        err = np.mean((pred - y) ** 2)
        if err < best_err:
            best_err = err
            best_init = np.array([a, alpha, beta, b], dtype=float)

    # fallback if grid search failed
    if best_init is None:
        best_init = np.array([1.0, 0.5, 0.2, 0.0], dtype=float)

    # refine all parameters with bounded L-BFGS-B on original MSE
    def objective(p):
        pred = scaling_law_func(X, p)
        return np.mean((pred - y) ** 2)

    # bounds: a>0, alpha>=0, beta>=0, 0<=b<=min(y)
    bounds = [
        (1e-12, None),    # a
        (0.0, None),      # alpha
        (0.0, None),      # beta
        (0.0, y_min)      # b
    ]

    res = minimize(objective, best_init, method='L-BFGS-B', bounds=bounds)
    if res.success:
        return res.x
    else:
        return best_init

Parallel Scaling Law

All Runs (sorted by R²)