Прогноз осадков с помощью машинного обучения

Сегодня метеорологи подсчитали, что 90% прогнозов погоды оказываются верными в течение пяти дней. Прогнозы обычно основаны на двух отдельных методах:

Подходы, основанные на физике. Эти подходы используют модели, которые интегрируют измеримые величины, такие как давление, движение облаков, состояние неба… Такие модели хороши для прогнозирования погоды на ближайшие дни или недели.
Подходы, не связанные с физикой (основанные на данных). Эти подходы используют исторические данные для создания моделей, которые могут делать прогнозы. Такие модели показывают хорошие результаты в прогнозировании погоды на срок до 6 часов или так называемом прогнозировании текущей погоды.

В этой статье мы обсудим вторую категорию подходов. Мы будем обсуждать различные форматы данных о погоде, как можно использовать машинное обучение для прогнозирования текущей погоды и как разработанные архитектуры могут быть полезны для решения аналогичных проблем.

Данные для прогнозирования климата

Поскольку подходы, свободные от физики, используют исторические данные, давайте начнем с изучения доступных данных.

Мы будем использовать два основных источника данных:

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

Табличные данные: эти данные имеют форму записей измеримых величин, таких как температура, влажность или скорость ветра.

Хотя оба источника данных важны для построения мощных моделей, мы сосредоточимся на первом (данные изображения, полученные с радаров или спутников) по соображениям простоты. Наиболее часто используемые модели с данными изображения: Сверточные нейронные сети (CNN. ).

После этой работы мы собираемся использовать архитектуру U-Net для создания нашей собственной модели прогнозирования текущей погоды.

Архитектура

Начинать с существующей архитектуры полезно по многим причинам, среди которых:

Архитектуры служат руководством для создания новых моделей. Люди, создающие новые архитектуры, применяют метод проб и ошибок, чтобы достичь окончательных результатов. Повторно используя их окончательные результаты, мы можем сэкономить много времени.

Предварительно обученные модели обычно доступны для немедленного использования. Когда исследователи публикуют свои новые архитектуры, они обычно также публикуют обученные параметры, чтобы пользователям не приходилось проходить через обучение/оптимизацию с нуля. Это особенно полезно для очень крупных, ресурсоемких моделей.

Примеры известных архитектур видения включают в себя:

LeNet (60 тыс. параметров)

AlexNet (60м параметров)

ВГГ-16 (параметры 138м)

U-Net

U-Net — это архитектура, основанная на полностью сверточной сети, то есть в ней нет полностью связанных слоев. Впервые он был представлен для задачи сегментации медицинских изображений. Эти результаты вдохновили исследователей распространить его на другие задачи компьютерного зрения.

В 2019 году Google использовал архитектуру на базе U-Net для создания модели прогнозирования осадков.

Название «U-Net» происходит от U-образной формы его архитектуры.

Мы выделяем три основных компонента:

Сжатие/кодирование: последовательность слоев свертки/объединения, которые сжимают входное изображение в представление меньшего размера.
Мост/узкое место: нижняя часть буквы «U» соединяет кодер с декодером. Он формируется серией операций свертки.
Деконтрактация/декодер: последовательность восходящих сверток и слоев свертки, которые «распаковывают» выходные данные узкого места.

Архитектура U-Net до сих пор выглядит как автокодировщик. Однако разница заключается в передаче информации между кодером и декодером. Эта информация передается путем объединения результатов свертки кодера с результатами свертки декодера. Эта модификация повышает разрешение и позволяет модели выводить более точные выходные данные (расположение пикселей в выходных данных будет более точным, что решает проблему, которая возникла в моделях, которые не включали конкатенацию). Объединение осуществляется симметрично.

Построение модели U-Net для прогнозирования текущих осадков

В этом разделе мы построим модель U-Net для прогнозирования погоды по спутниковым изображениям. Мы будем использовать предварительно обученную модель под названием Rain-Net.

Код ниже доступен в этой colab.

Сначала мы устанавливаем wradlib, библиотеку с открытым исходным кодом для обработки данных метеорадиолокаторов.

!pip install wradlib
import wradlib as wrl

Затем мы пишем две служебные функции для загрузки спутниковых данных с сервера открытых данных DWD.


import urllib.request
import io
import numpy as np
import datetime

def download_and_read_RY(RY_timestamp):
   url = f"https://opendata.dwd.de/weather/radar/radolan/ry/raa01-ry_10000-{RY_timestamp}-dwd---bin"
   data_binary = urllib.request.urlopen(url).read()
   data, attr = wrl.io.read_radolan_composite(io.BytesIO(data_binary), missing = 0)
   data = data.astype("float32")
   return data,attr

def download_data():
   latest_scan, latest_attr = download_and_read_RY("latest")
   latest_datetime = latest_attr["datetime"]
   list_for_downloading = [ts.strftime("%y%m%d%H%M") for ts in
       [latest_datetime - datetime.timedelta(minutes=t) for t in [15, 10, 5]]]
   previous_scans = np.array([download_and_read_RY(ts)[0] for ts in list_for_downloading])
   print(list_for_downloading)
   scans = np.concatenate([previous_scans, latest_scan[np.newaxis, ::, ::]], axis = 0)

   return scans, latest_datetime

Эти утилиты позволяют нам загружать последние 4 спутниковых изображения — именно столько спутниковых изображений нам нужно для прогнозирования с помощью нашей предварительно обученной модели.

Затем мы можем использовать созданные функции, чтобы получить последние 4 изображения.


RY_latest, RY_latest_timestep = download_data()

После получения изображений мы используем метод wradlib vis.plot_ppi для построения данных.

for i in range(RY_latest.shape[0]):
   wrl.vis.plot_ppi(RY_latest[i])

Теперь мы загрузили наши данные. Давайте загрузим модель дальше.

Начнем с импорта соответствующих классов. В этой статье мы будем использовать TensorFlow.

из tensorflow.keras.layers import Input, Conv2D, Activation, Concatenate, Conv2DTranspose, MaxPool2D

из tensorflow.keras.models Модель импорта

Давайте построим 3 примитивных строительных блока. Эти «строительные блоки» будут использоваться для создания всей архитектуры согласно этой реализации.

Первый блок соответствует последовательности сверточных слоев, мы называем его «conv_block».

def conv_block(input, num_filters):
   x = Conv2D(num_filters, 3, padding ="same", kernel_initializer="he_normal"(input))
   x = Activation("relu")(x)
   x = Conv2D(num_filters, 3, padding ="same", kernel_initializer="he_normal"(x))
   x = Activation("relu")(x)
   return x

Второй блок используется для построения части кодера (сверточный блок + максимальный пул). Мы называем это «encoder_block».

def encoder_block(input, num_filters):
   x = conv_block(input, num_filters)
   p = MaxPool2D((2,2))(x)
   return x,p

Третий и последний блок — это «decoder_block» (свертка + конкатенация + свертка).

def decoder_block(input, skip_features, num_filters):
   x = Conv2DTranspose(num_filters, (2,2), strides=2, padding="same", kernel_initializer="he_normal")
   x = Concatenate()([x, skip_features])
   x = conv_block(x, num_filters)
   return x

Мы объединяем эти строительные блоки, чтобы построить модель U-Net.

def build_unet(input_shape = (928, 928, 4)):

   inputs = Input(input_shape)
   e1, p1 = encoder_bock(inputs, 64)
   e2, p2 = encoder_bock(p1, 128)
   e3, p3 = encoder_bock(p2, 256)
   e4, p4 = encoder_bock(p3, 512)

   b1 = conv_block(p4, 1024)

   d1 = decoder_block(b1, e4, 512)
   d2 = decoder_block(d1, e3, 256)
   d3 = decoder_block(d2, e2, 128)
   d4 = decoder_block(d3, e1, 64)

   outputs = Conv2D(1, 1, padding="same", activation="linear")(d4)

   unet = Model(inputs, outputs, name="U-Net")

   return unet

Не стесняйтесь экспериментировать с реализацией, улучшать ее или адаптировать под свои нужды, если вы хотите использовать ее для чего-то еще.

Обучение этой модели может занять много времени. К счастью, существует модель Rain-Net, созданная на основе архитектуры U-Net и специализирующаяся на прогнозировании осадков.

Давайте клонируем его репозиторий GitHub.

! git clone https://github.com/hydrogo/rainnet.git

Затем мы загружаем предварительно обученные веса для этой модели.

!wget -O /content/rainnet/rainnet_weights.h5 https://zenodo.org/record/3630429/files/rainnet_weights.h5

Следующим шагом будет создание модели на основе архитектуры, найденной в репозитории, а затем загрузка весов, загруженных в эту модель.

import sys
from rainnet import rainnet
model = rainnet.rainnet()
model.load_weights('/content/rainnet/rainnet_weights.h5')
model.summary()

Загруженные нами изображения имеют размер 900*900 пикселей. Мы собираемся изменить форму этих изображений, чтобы они соответствовали ожидаемым данным Rain-Net.

def Scaler(array):
   return np.log(array+0.01)


def invScaler(array):
   return np.exp(array) - 0.01


def pad_to_shape(array, from_shape=900, to_shape=928, how="mirror"):
   # calculate how much to pad in respect with native resolution
   padding = int( (to_shape - from_shape) / 2)
   # for input shape as (batch, W, H, channels)
   if how == "zero":
       array_padded = np.pad(array, ((0,0),(padding,padding),(padding,padding),(0,0)), mode="constant", constant_values=0)
   elif how == "mirror":
       array_padded = np.pad(array, ((0,0),(padding,padding),(padding,padding),(0,0)), mode="reflect")
   return array_padded


def pred_to_rad(pred, from_shape=928, to_shape=900):
   # pred shape 12,928,928
   padding = int( (from_shape - to_shape) / 2)
   return pred[::, padding:padding+to_shape, padding:padding+to_shape].copy()

'''
the spatial extent of input data has to be a multiple of 2n+1
where n is the number of max pooling layers
'''

def data_preprocessing(X):

   # 0. Right shape for batch
   X = np.moveaxis(X, 0, -1)
   X = X[np.newaxis, ::, ::, ::]
   # 1. To log scale
   '''
   We use a log scale to respond to skewness towards large values
   '''
   X = Scaler(X)
   # 2. from 900x900 to 928x928
   X = pad_to_shape(X)

   return X


def data_postprocessing(nwcst):

   # 0. Squeeze empty dimensions
   nwcst = np.squeeze(np.array(nwcst))
   # 1. Convert back to rainfall depth
   nwcst = invScaler(nwcst)
   # 2. Convert from 928x928 back to 900x900
   nwcst = pred_to_rad(nwcst)
   # 3. Return only positive values
   nwcst = np.where(nwcst>0, nwcst, 0)
   return nwcst

Затем мы создаем функцию, которая делает прогнозы.

def prediction(model_instance, input_data, lead_time=24):

   input_data = data_preprocessing(input_data)
   nwcst = []

   for _ in range(lead_time):
       pred = model_instance.predict(input_data)
       nwcst.append(pred)
       input_data = np.concatenate([input_data[::, ::, ::, 1:], pred], axis = -1)

   nwcst = data_postprocessing(nwcst)

   return nwcst

Затем мы вызываем эту функцию для данных, которые мы скачали ранее.

Y_pred = prediction(model, RY_latest)

Мы можем построить прогнозы и сохранить их, чтобы использовать сохраненные результаты для создания изображения в формате gif, которое позволяет нам визуализировать прогнозы.

import matplotlib.pyplot as plt
names = []
for i in range(0, 19):
   im = wrl.vis.plot_ppi(Y_pred[i])
   name = '/content/image'+str(i)+'.png'
   names.append(name)
   plt.savefig(name)

import imageio
from google.colab import files

imgs = []

for name in names:
   imgs.append(imageio.imread(name))

imageio.mimsave('/content/pred.gif', imgs, fps=4)
files.download('/content/pred.gif')

Поздравляем с тем, что вы зашли так далеко! Теперь вы можете использовать Rain-Net, чтобы делать прогнозы и визуализировать их.

Заключение

В этой статье мы использовали модель машинного обучения (Rain-Net) для прогнозирования текущей погоды. Мы использовали архитектуру U-Net, которую построили с помощью TensorFlow. Мы загрузили предварительно обученную модель для прогнозирования спутниковых изображений.

Эту реализацию можно улучшить разными способами. Например:

Точная настройка модели в вашем наборе данных.
Используйте в архитектуре модуль, основанный на внимании, например CBAM (Модуль сверточного блока).

Использованная литература

Приходите на один из наших бесплатных мастер-классов!

Начните свою карьеру в качестве специалиста по данным с наших бесплатных семинаров, которые основаны на адаптируемой учебной программе и проводятся под руководством отраслевых экспертов.