Примечание
Перейти к концу, чтобы загрузить полный пример кода. или запустить этот пример в вашем браузере через JupyterLite или Binder
Прогнозирование уровня CO2 по набору данных Mona Loa с помощью регрессии Гауссова процесса (GPR)
Этот пример основан на разделе 5.4.3 «Гауссовы процессы для машинного обучения» [1]. Он демонстрирует пример сложного проектирования ядра и оптимизации гиперпараметров с использованием градиентного восхождения по логарифмическому маргинальному правдоподобию. Данные состоят из ежемесячных средних концентраций атмосферного CO2 (в частях на миллион по объему (ppm)), собранных в обсерватории Мауна-Лоа на Гавайях, в период с 1958 по 2001 год. Цель состоит в том, чтобы смоделировать концентрацию CO2 как функцию времени \(t\) и экстраполировать на годы после 2001 года.
Ссылки
print(__doc__) # Authors: The scikit-learn developers # SPDX-License-Identifier: BSD-3-Clause
Создание набора данных
Мы получим набор данных из обсерватории Мауна-Лоа, которая собирала образцы воздуха. Мы заинтересованы в оценке концентрации CO2 и экстраполяции ее на будущие годы. Сначала мы загружаем исходный набор данных, доступный в OpenML, в виде pandas DataFrame. Это будет заменено Polars, как только fetch_openml добавит в него родную поддержку.
from sklearn.datasets import fetch_openml co2 = fetch_openml(data_id=41187, as_frame=True) co2.frame.head()
Сначала мы обрабатываем исходный DataFrame, чтобы создать столбец даты и выбрать его вместе со столбцом CO2.
import polars as pl
co2_data = pl.DataFrame(co2.frame[["year", "month", "day", "co2"]]).select(
pl.date("year", "month", "day"), "co2"
)
co2_data.head()
co2_data["date"].min(), co2_data["date"].max()
(datetime.date(1958, 3, 29), datetime.date(2001, 12, 29))
Мы видим, что получаем концентрацию CO2 в некоторые дни с марта 1958 года по декабрь 2001 года. Мы можем построить эту исходную информацию, чтобы лучше понять ее.
import matplotlib.pyplot as plt
plt.plot(co2_data["date"], co2_data["co2"])
plt.xlabel("date")
plt.ylabel("CO$_2$ concentration (ppm)")
_ = plt.title("Raw air samples measurements from the Mauna Loa Observatory")
Мы предобработаем набор данных, взяв ежемесячное среднее значение и удалив месяцы, для которых не были собраны измерения. Такая обработка сгладит данные.
co2_data = (
co2_data.sort(by="date")
.group_by_dynamic("date", every="1mo")
.agg(pl.col("co2").mean())
.drop_nulls()
)
plt.plot(co2_data["date"], co2_data["co2"])
plt.xlabel("date")
plt.ylabel("Monthly average of CO$_2$ concentration (ppm)")
_ = plt.title(
"Monthly average of air samples measurements\nfrom the Mauna Loa Observatory"
)
В этом примере идея заключается в прогнозировании концентрации CO2 в зависимости от даты. Мы также заинтересованы в экстраполяции на последующие годы после 2001 года.
В качестве первого шага мы разделим данные и целевую переменную для оценки. Поскольку данные представляют собой дату, мы преобразуем ее в числовое значение.
X = co2_data.select(
pl.col("date").dt.year() + pl.col("date").dt.month() / 12
).to_numpy()
y = co2_data["co2"].to_numpy()
Проектирование подходящего ядра
Для проектирования ядра, используемого с нашим гауссовым процессом, мы можем сделать некоторые предположения относительно имеющихся данных. Мы наблюдаем, что у них есть несколько характеристик: мы видим долгосрочную возрастающую тенденцию, выраженную сезонную вариацию и некоторые незначительные нарушения. Мы можем использовать различные подходящие ядра, которые отражали бы эти особенности.
Во-первых, долгосрочную тенденцию роста можно подобрать с помощью ядра радиальной базисной функции (RBF) с большим параметром длины масштаба. Ядро RBF с большим масштабом длины обеспечивает плавность этого компонента. Возрастающая тенденция не навязывается, чтобы предоставить свободу нашему моделированию. Конкретный масштаб длины и амплитуда являются свободными гиперпараметрами.
from sklearn.gaussian_process.kernels import RBF long_term_trend_kernel = 50.0**2 * RBF(length_scale=50.0)
Сезонная вариация объясняется периодическим экспоненциальным синус-квадратным ядром с фиксированной периодичностью 1 год. Масштаб длины этого периодического компонента, контролирующий его плавность, является свободным параметром. Для того, чтобы позволить затухание при отклонении от точной периодичности, берется произведение с ядром RBF. Масштаб длины этого компонента RBF управляет временем затухания и является дополнительным свободным параметром. Такой тип ядра также известен как локально периодическое ядро.
from sklearn.gaussian_process.kernels import ExpSineSquared
seasonal_kernel = (
2.0**2
* RBF(length_scale=100.0)
* ExpSineSquared(length_scale=1.0, periodicity=1.0, periodicity_bounds="fixed")
)
Незначительные нарушения необходимо объяснить компонентом ядра рациональной квадратичной функции, масштаб длины и параметр альфа которого, количественно характеризующий диффузность масштабов длины, необходимо определить. Ядро рациональной квадратичной функции эквивалентно ядру RBF с несколькими масштабами длины и лучше учитывает различные нарушения.
from sklearn.gaussian_process.kernels import RationalQuadratic irregularities_kernel = 0.5**2 * RationalQuadratic(length_scale=1.0, alpha=1.0)
Наконец, шум в наборе данных можно учесть с помощью ядра, состоящего из вклада ядра RBF, который должен объяснить коррелированные составляющие шума, такие как местные погодные явления, и вклада белого ядра для белого шума. Относительные амплитуды и масштаб длины RBF являются дополнительными свободными параметрами.
from sklearn.gaussian_process.kernels import WhiteKernel
noise_kernel = 0.1**2 * RBF(length_scale=0.1) + WhiteKernel(
noise_level=0.1**2, noise_level_bounds=(1e-5, 1e5)
)
Таким образом, наше окончательное ядро является суммой всех предыдущих ядер.
co2_kernel = (
long_term_trend_kernel + seasonal_kernel + irregularities_kernel + noise_kernel
)
co2_kernel
50**2 * RBF(length_scale=50) + 2**2 * RBF(length_scale=100) * ExpSineSquared(length_scale=1, periodicity=1) + 0.5**2 * RationalQuadratic(alpha=1, length_scale=1) + 0.1**2 * RBF(length_scale=0.1) + WhiteKernel(noise_level=0.01)
Подгонка модели и экстраполяция
Теперь мы готовы использовать регрессор гауссова процесса и подогнать имеющиеся данные. Чтобы следовать примеру из литературы, мы вычтем среднее значение из целевой переменной. Мы могли бы использовать normalize_y=True. Однако в этом случае мы также масштабировали бы целевую переменную (разделив y на ее стандартное отклонение). Таким образом, гиперпараметры различных ядер имели бы разное значение, поскольку они не были бы выражены в ppm.
from sklearn.gaussian_process import GaussianProcessRegressor y_mean = y.mean() gaussian_process = GaussianProcessRegressor(kernel=co2_kernel, normalize_y=False) gaussian_process.fit(X, y - y_mean)
Теперь мы воспользуемся гауссовым процессом для прогнозирования:
- данных обучения для проверки качества подгонки;
- будущих данных для просмотра экстраполяции, выполняемой моделью.
Таким образом, мы создаем синтетические данные с 1958 года по текущий месяц. Кроме того, нам нужно добавить вычтенное среднее значение, вычисленное во время обучения.
import datetime import numpy as np today = datetime.datetime.now() current_month = today.year + today.month / 12 X_test = np.linspace(start=1958, stop=current_month, num=1_000).reshape(-1, 1) mean_y_pred, std_y_pred = gaussian_process.predict(X_test, return_std=True) mean_y_pred += y_mean
plt.plot(X, y, color="black", linestyle="dashed", label="Measurements")
plt.plot(X_test, mean_y_pred, color="tab:blue", alpha=0.4, label="Gaussian process")
plt.fill_between(
X_test.ravel(),
mean_y_pred - std_y_pred,
mean_y_pred + std_y_pred,
color="tab:blue",
alpha=0.2,
)
plt.legend()
plt.xlabel("Year")
plt.ylabel("Monthly average of CO$_2$ concentration (ppm)")
_ = plt.title(
"Monthly average of air samples measurements\nfrom the Mauna Loa Observatory"
)
Наша обученная модель способна правильно подогнать предыдущие данные и экстраполировать на будущие годы с высокой уверенностью.
Интерпретация гиперпараметров ядра
Теперь мы можем взглянуть на гиперпараметры ядра.
gaussian_process.kernel_
44.8**2 * RBF(length_scale=51.6) + 2.64**2 * RBF(length_scale=91.5) * ExpSineSquared(length_scale=1.48, periodicity=1) + 0.536**2 * RationalQuadratic(alpha=2.89, length_scale=0.968) + 0.188**2 * RBF(length_scale=0.122) + WhiteKernel(noise_level=0.0367)
Таким образом, большая часть целевого сигнала, с вычтенной средней величиной, объясняется долгосрочной восходящей тенденцией на ~45 ppm и длиной масштаба ~52 года. Периодический компонент имеет амплитуду ~2,6 ppm, время затухания ~90 лет и длину масштаба ~1,5. Большое время затухания указывает на то, что у нас есть компонент, очень близкий к сезонной периодичности. Коррелированный шум имеет амплитуду ~0,2 ppm с длиной масштаба ~0,12 года и вклад белого шума ~0,04 ppm. Таким образом, общий уровень шума очень невелик, что указывает на то, что данные могут быть очень хорошо объяснены моделью.
Общее время выполнения скрипта: (0 минут 4,187 секунды)
Связанные примеры
© 2007–2025 The scikit-learn developers
Licensed under the 3-clause BSD License.
https://scikit-learn.org/1.6/auto_examples/gaussian_process/plot_gpr_co2.html
