Problem 445070 · medium · Level 04 Non-Linear Data Structures

Which Peak Is the Right One?

non-convex loss · local minima · random restarts · gradient descent · chain rule

A spectrometer shows several bell-shaped peaks, and you want the position c of the strongest one. You fit a single peak of height 1 and known width width,

m_c(x) = exp(-(x - c)² / (2 · width²))

by minimising the mean squared error L(c) = mean over the data of (m_c(x) - y)² over the one parameter c. The loss has a valley near every peak, and far from the data it is almost flat.

Write find_peak(xs, ys, width, restarts, lr, steps, seed):

  1. Create rng = random.Random(seed) and draw the starting points one after another, each as rng.uniform(min(xs), max(xs)), restarts of them.
  2. From every start, run steps steps of gradient descent c ← c - lr · L'(c).
  3. Return (best_c, best_loss, finals): finals is the list of final positions in the order of the starts, and best_c is the final position with the lowest loss, best_loss its loss.

The helper spectrum(peaks, width, lo, hi, k) builds test data: k equally spaced positions from lo to hi, with the heights of the peaks (centre, height) added up and rounded to 3 decimals.

Examples

Input:  xs, ys = spectrum([(3, 1.0), (7.5, 0.6)], 0.6, 0, 10, 41)
        find_peak(xs, ys, 0.6, 6, 0.5, 300, 0)
Output: (3.000001, 0.037355, [7.499997, 7.499997, 3.000001, 3.000001, 3.000001, 3.000001])   (rounded here)
Explanation: two starts (8.44 and 7.58) roll into the valley of the weaker peak at 7.5,
the other four find the strong peak at 3. The best loss is at c = 3.

Input:  the same data, find_peak(xs, ys, 0.6, 6, 0.5, 300, 2)
Output: (3.000001, 0.037355, [11.366716, 11.35867, -1.363505, 3.000001, 7.499997, 7.499997])   (rounded here)
Explanation: starts near the ends of the range (9.56, 9.48, 0.57) drift away from all the data:
there the loss only falls by moving the model peak further out.

Constraints

  • 5 <= len(xs) == len(ys) <= 60, 1 <= restarts <= 10, 1 <= steps <= 300, 0 < lr <= 1, 0.3 <= width <= 2
  • draw exactly one rng.uniform(min(xs), max(xs)) per restart, all from the same rng, in order
  • floats are compared with a tolerance of 1e-6; the tests have one clearly best final position

Goals

  • Derive the gradient of a loss whose model is not linear in its parameter
  • See gradient descent settle in different valleys, or drift on a flat plateau, depending on the start
  • Use seeded random restarts and keep the best result
Starting Python…