Подгонка логистической функции к ВВП Китая: почему высокое R² не спасает прогноз 1970

Разбираем подгонку логистической кривой к ВВП Китая с scipy.optimize.curve_fit: почему нормировка и высокий R² обманывают и прогноз 1970 даёт большую абсолютную ошибку.

На первый взгляд логистическая кривая, подогнанная к макроэкономическим временным рядам, выглядит убедительно: функция гладкая, R² высокий, и наложенная кривая хорошо следует данным. Однако стоит сделать прогноз на один конкретный год и вернуть масштаб к исходным единицам, как абсолютная ошибка может оказаться большой. Ниже приведено практическое объяснение, почему так происходит, что именно делает код, и несколько безопасных способов рассуждать об этом без лишних допущений.

Постановка задачи

Нужно подогнать логистическую функцию к ряду ВВП Китая с помощью scipy.optimize.curve_fit, предварительно нормировав признаки и целевые значения к [0, 1]. Кривая визуально совпадает с данными, а отчетный R² очень высок. Но после обратного преобразования масштабов прогноз для 1970 года получается около 1.92285528141e11, тогда как фактический ВВП за тот год — примерно 9.15062113063745e10.

Воспроизводимый пример кода

Ниже скрипт, который нормирует входы и выходы, подбирает сигмоиду, считает R² на отложенном срезе и получает прогноз для 1970 года после обратного масштабирования. Имена переменных выбраны ради ясности, логика соответствует описанному сценарию.

import numpy as np
import pandas as pd
from scipy.optimize import curve_fit
from scipy.special import expit

frame = pd.read_csv("china_gdp.csv")
split_mask = np.random.rand(len(frame)) < 0.8

t_raw = frame['Year'].values.astype(float)
g_raw = frame['Value'].values.astype(float)

# масштабирование к [0, 1]
t_unit = (t_raw - t_raw.min()) / (t_raw.max() - t_raw.min())
g_unit = (g_raw - g_raw.min()) / (g_raw.max() - g_raw.min())

def sig_curve(x, a1, a2, a3, a4):
    return a3 + a4 * expit(a1 * (x - a2))

init_guess = (5.0, 0.5, 0.0, 1.0)
param_bounds = ([0.0, 0.0, -np.inf, 0.0], [np.inf, 1.0, np.inf, np.inf])

opt_params, covar = curve_fit(sig_curve, t_unit, g_unit, p0=init_guess, bounds=param_bounds)
print(dict(zip(['a1','a2','a3','a4'], opt_params)))

import matplotlib.pyplot as plt

t_grid = np.linspace(0, 1, 300)
g_grid = sig_curve(t_grid, *opt_params)
plt.scatter(t_unit, g_unit, s=15, label='data')
plt.plot(t_grid, g_grid, label='fit')
plt.legend(); plt.show()

from sklearn.metrics import r2_score

# примечание: подгонка выше использовала все точки; этот срез нужен лишь для иллюстрации
thold = t_unit[~split_mask]
ghold = g_unit[~split_mask]
gpred = sig_curve(thold, *opt_params)
print(r2_score(ghold, gpred))

yr_query = 1970
scaled_pred = sig_curve((yr_query - t_raw.min()) / (t_raw.max() - t_raw.min()), *opt_params)
back_to_units = scaled_pred * (g_raw.max() - g_raw.min()) + g_raw.min()
print(back_to_units)

Типичные результаты для этой настройки: R² порядка 0.9985891574981185 на нормированном отложенном срезе и прогноз для 1970 года около 192285528141.3661, тогда как фактическое значение — 91506211306.3745.

Почему на вид отличная подгонка промахивается по конкретному году

С реализацией всё в порядке. Расхождение — следствие взаимодействия выбранного класса функций с данными. Модель представляет собой «константа плюс сигмоида». Сигмоида монотонна, имеет пологие края и крутой средний участок. В рассматриваемом ряду ВВП есть выраженная фаза ускорения в начале–середине 1990-х. Когда и годы, и значения ВВП нормируются в [0, 1], оптимизация идет в масштабированном пространстве, а не в исходных единицах. Это делает относительные ошибки в нижней части кривой небольшими в нормированных терминах, хотя после обратного преобразования в доллары они могут оборачиваться крупными абсолютными отклонениями. Одной сигмоиды просто недостаточно, чтобы точно описать все изгибы и сдвиги кривой, где присутствует резкий период ускорения.

Отсюда и то, почему суммарная ошибка может выглядеть умеренной в шкале [0, 1]. Например, сумма модулей ошибок на нормированных выходах может быть около 0.26858. Квадратичные метрики в [0, 1] «сжимают» большие различия в исходных единицах, поэтому они недооценивают абсолютные отклонения в нижней части диапазона после обратного масштабирования.

Есть и еще одна тонкость. Вызов подгонки выше использует все точки для оценки параметров, а затем метрика считается на подмножестве, взятом из той же временной серии. Даже если переобучать только на замаскированной части, случайное разбиение временного ряда — слабая стратегия валидации для прогноза. В любом случае корневая проблема здесь в том, что класс функций слишком жесткий, чтобы одновременно «обнять» и низкие значения ВВП, и участок бурного роста после нормировки.

Практичный способ уточнить прогнозы ранних лет

Ограничение подгонки режимом, где нет резкого взлета, согласует форму функции с локальным поведением данных. Усечение ряда до 1960–1991 годов в данном случае заметно улучшает оценку для 1970 года. Ниже показано, как выполнить подгонку на этом раннем окне.

_take = 32
split_mask = np.random.rand(_take) < 0.8

t_slice = frame['Year'].values.astype(float)[:_take]
g_slice = frame['Value'].values.astype(float)[:_take]

# нормировка внутри среза
ts = (t_slice - t_slice.min()) / (t_slice.max() - t_slice.min())
gs = (g_slice - g_slice.min()) / (g_slice.max() - g_slice.min())

opt_params2, covar2 = curve_fit(sig_curve, ts, gs, p0=init_guess, bounds=param_bounds)

yr_query = 1970
scaled_q = (yr_query - t_slice.min()) / (t_slice.max() - t_slice.min())
scaled_ans = sig_curve(scaled_q, *opt_params2)
back_scaled = scaled_ans * (g_slice.max() - g_slice.min()) + g_slice.min()
print(back_scaled)

Оценка для 1970 года по усеченной подгонке получается около 94736151945.78181 — гораздо ближе к 91506211306.3745. Визуально кривая на раннем сегменте также плотнее прилегает к области низких значений.

Как исследовать гибкость модели, не скатываясь сразу в переобучение

Если цель — системно посмотреть, как добавление свободы меняет качество подгонки, можно поэкспериментировать с полиномиальной базовой моделью, варьируя степень, и наблюдать, как прогнозы переходят от недообучения к переобучению по мере добавления параметров. Фрагмент ниже выполняет МНК-подгонку для нескольких степеней полинома и строит прогнозы относительно фиксированного «реального» вектора. Это не рецепт моделирования ВВП; это компактный способ нащупать интуицию о числе степеней свободы и их влиянии на поведение подгонки.

import numpy as np
import matplotlib.pyplot as plt

step_deg = 3
rng_seed = 1

orders = step_deg * (1 + np.arange(6))
count = 50
x_ax = np.arange(count)
np.random.seed(int(rng_seed)); y_true = 1e2 * np.random.rand(count)

y_curves = [f(x_ax) for f in [np.poly1d(np.polyfit(x_ax, y_true, k)) for k in orders]]

fig, axs = plt.subplots(3, 2)
slots = sorted((i % 3, i % 2) for i in range(6))
for i in range(len(orders)):
    axs[*slots[i]].set_title(f'Degree = {orders[i]}')
    axs[*slots[i]].scatter(x_ax, y_true, s=15, label='data')
    axs[*slots[i]].plot(x_ax, y_curves[i], label='fit')
    axs[*slots[i]].legend()
plt.show()

Это упражнение наглядно показывает, как дополнительные параметры сначала «подтягивают» подгонку, а затем, перейдя некий порог, начинают следовать каждой «рифлености» данных. Оно также подчеркивает, почему одной сигмоиды может быть мало для данных, в которых спокойные периоды чередуются с взрывным ростом.

Почему это важно

Выбор формы функции, соответствующей природе явления, не менее важен, чем оптимизатор и метрика. Нормировка целевых переменных в [0, 1] — привычная и полезная практика, но она переносит оптимизацию в пространство, где доминируют относительные различия. После обратного масштабирования кажущиеся малыми расхождения в нормированной шкале могут превратиться в крупные абсолютные ошибки в нижнем диапазоне. Если затем интерпретировать их в исходных единицах без учета этого контекста, выводы могут быть обманчивыми. Важен и протокол валидации: случайные разбиения во временных рядах не отражают реальное вневыборочное прогнозирование, а подгонка на всех точках с последующим отчетом по разбиению может завышать оценку качества.

Итоги

Когда модель дает отличный общий фит, но промахивается по отдельным годам в абсолютных величинах, сначала проверьте согласование класса функций с режимами данных. Для ВВП с резким взлетом одна логистическая кривая распределяет ошибки по всему нормированному диапазону и не сможет одновременно точно попасть и в ранние низкие значения, и в поздние высокие. Подгонка в согласованном режиме, как в срезе 1960–1991, может заметно улучшить прогнозы ранних лет. Чтобы понять, как емкость влияет на качество подгонки, варьируйте число степеней свободы в простой базовой модели и наблюдайте, где проходит граница между здоровой гибкостью и переобучением. И главное — интерпретируйте ошибки в той шкале, которая для вас важна, и осторожно относитесь к случайным разбиениям временных рядов при оценке прогностической надежности.