Precipitation Nowcasting so strojovým učením

V súčasnosti meteorológovia odhadujú, že 90% predpovedí počasia je správnych v priebehu 5 dní. Predpovede, ktoré sa robia, sú zvyčajne založené na dvoch samostatných metódach:

Fyzikálne prístupy: Tieto prístupy využívajú modely, ktoré integrujú merateľné veličiny ako tlak, pohyb mrakov, podmienky na oblohe... Takéto modely sú dobré na predpovedanie počasia na najbližšie dni alebo týždne.
Prístupy bez fyziky (založené na údajoch): Tieto prístupy využívajú historické údaje na vytváranie modelov, ktoré môžu predpovedať. Takéto modely vykazujú dobré výsledky pri predpovedaní počasia až na 6 hodín alebo čo je známe ako weather nowcasting.

V tomto článku sa budeme zaoberať druhou kategóriou prístupov. Budeme diskutovať o rôznych formátoch údajov o počasí, o tom, ako možno využiť strojové učenie pri vykonávaní aktuálneho počasia a ako môžu byť vyvinuté architektúry prospešné pri riešení podobných problémov.

Údaje pre predpoveď klímy

Keďže prístupy bez fyziky používajú historické údaje, začnime pohľadom na dostupné údaje.

Budeme používať dva hlavné zdroje údajov:

Obrazové údaje: Tieto údaje majú formu radarových alebo satelitných snímok špecifikovanej geografickej oblasti. Používa sa na predpovedanie zrážok, pohybu vetra alebo vlhkosti.

Tabuľkové údaje : Tieto údaje majú formu záznamov merateľných veličín ako je teplota, vlhkosť alebo rýchlosť vetra.

Aj keď sú oba zdroje údajov dôležité na vytvorenie výkonných modelov, z dôvodov jednoduchosti sa zameriame na prvé (obrazové údaje zozbierané z radarov alebo satelitov). Najbežnejšie používané modely s obrazovými údajmi sú Convolutional Neural Networks (CNN ).

Po tejto práci použijeme architektúru U-Net na vytvorenie vlastného modelu nowcastingu.

Architektúra

Vychádzať z existujúcej architektúry je užitočné z mnohých dôvodov, medzi ktoré patrí:

Architektúry slúžia ako návod na vytváranie nových modelov. Ľudia, ktorí vytvárajú nové architektúry, používajú metódu pokus-omyl, aby sa dostali ku konečným výsledkom. Opätovným použitím ich konečných výsledkov môžeme ušetriť veľa času.

Predškolené modely sú zvyčajne k dispozícii na okamžité použitie. Keď výskumníci zverejnia svoje nové architektúry, zvyčajne zverejnia aj natrénované parametre, takže používatelia nebudú musieť prechádzať problémami s tréningom/optimalizáciou od začiatku. Toto je obzvlášť užitočné pre veľmi veľké modely, ktoré si vyžadujú zdroje.

Príklady známych architektúr vízií zahŕňajú:

LeNet (60 000 parametrov)

AlexNet (60 m parametrov)

VGG-16 (138 m parametrov)

U-Net

U-Net je architektúra založená na plne konvolučnej sieti, čo znamená, že nemá žiadne úplne prepojené vrstvy. Prvýkrát bol predstavený pre úlohu segmentácie medicínskeho obrazu. Tieto výsledky inšpirovali výskumníkov k rozšíreniu na ďalšie úlohy v oblasti počítačového videnia.

V roku 2019 použil Google na vytvorenie modelu predpovede zrážok architektúru založenú na U-Net.

Názov „U-Net“ pochádza z tvaru „U“ jeho architektúry.

Rozlišujeme tri hlavné zložky:

Zmršťovanie / kódovanie: Postupnosť konvolučných/združovacích vrstiev, ktoré komprimujú vstupný obraz do reprezentácie menšej veľkosti.
Most/úzky profil: Spodná časť písmena „U“ spája kódovač s dekodérom. Tvorí ho séria konvolučných operácií.
Dekontrakcia / Dekodér : Postupnosť vzostupných a konvolučných vrstiev, ktoré „dekomprimujú“ výstup úzkeho miesta.

Architektúra U-Net zatiaľ vyzerá ako automatický kódovač. Rozdiel však spočíva v prenose informácií medzi kodérom a dekodérom. Tieto informácie sa odovzdávajú zreťazením výsledkov konvolúcií z kodéra s výsledkami konvolúcie dekodéra. Táto úprava zvyšuje rozlíšenie a umožňuje modelu vydávať priestorovo presnejší výstup (umiestnenie pixelov vo výstupe bude presnejšie, čím sa vyrieši problém, ktorý sa vyskytol v modeloch, ktoré nezahŕňali zreťazenie). Reťazenie sa robí symetricky.

Vytvorenie modelu U-Net pre zrážkové teraz

V tejto časti vytvoríme model U-Net na predpovedanie počasia zo satelitných snímok. Použijeme predtrénovaný model s názvom Rain-Net.

Nižšie uvedený kód je k dispozícii v tomto colab.

Najprv nainštalujeme wradlib, knižnicu s otvoreným zdrojom na spracovanie údajov meteorologických radarov

!pip install wradlib
import wradlib as wrl

Potom napíšeme dve pomocné funkcie na stiahnutie satelitných údajov z otvoreného dátového servera 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

Tieto nástroje nám umožňujú stiahnuť najnovšie 4 satelitné snímky, čo je počet satelitných snímok, ktoré potrebujeme na predpovedanie nášho vopred pripraveného modelu.

Potom môžeme použiť vytvorené funkcie na získanie najnovších 4 obrázkov


RY_latest, RY_latest_timestep = download_data()

Po získaní obrázkov použijeme metódu wradlib vis.plot_ppi na vykreslenie údajov

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

Teraz sme načítali naše údaje. Ďalej načítajme model.

Začneme importovaním príslušných tried. V tomto článku budeme používať TensorFlow.

z tensorflow.keras.layers importovať Input, Conv2D, Activation, Concatenate, Conv2DTransose, MaxPool2D

z tensorflow.keras.models importovať Model

Postavme 3 primitívne stavebné bloky. Tieto „stavebné bloky“ sa použijú na vytvorenie celej architektúry podľa tejto implementácie.

Prvý blok zodpovedá postupnosti konvolučných vrstiev, nazývame ho „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

Druhý blok sa používa na zostavenie časti kódovača (konvolučný blok + max pooling). Nazývame to „blok_kódovania“

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

Tretí a posledný blok je „decoder_block“ (upconvolúcia + zreťazenie + konvolúcia).

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

Tieto stavebné bloky kombinujeme na vytvorenie modelu 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

Pokojne sa pohrajte s implementáciou, vylepšite ju alebo si ju prispôsobte svojim potrebám, ak ju chcete použiť na niečo iné.

Tréning tohto modelu môže zabrať veľa času. Našťastie existuje model s názvom Rain-Net, ktorý bol vytvorený na základe architektúry U-Net a je špecializovaný na zrážkové odlievanie.

Poďme naklonovať jeho úložisko GitHub

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

Potom stiahneme predtrénované závažia pre tento model

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

Ďalším krokom je vytvorenie modelu založeného na architektúre nájdenej v úložisku a potom načítanie váh stiahnutých do tohto modelu

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

Obrázky, ktoré sme stiahli, majú veľkosť 900*900 pixelov. Tieto obrázky pretvoríme tak, aby zodpovedali očakávanému vstupu 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

Potom vytvoríme funkciu, ktorá robí predpovede.

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

Túto funkciu potom voláme na dátach, ktoré sme si predtým stiahli

Y_pred = prediction(model, RY_latest)

Predpovede môžeme vykresliť a uložiť, aby sme použili uložené výsledky na vytvorenie obrázka gif, ktorý nám umožní vizualizovať predpovede

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')

Gratulujeme, že ste to dotiahli až sem! Teraz môžete použiť Rain-Net na vytváranie predpovedí a ich vizualizáciu.

Záver

V tomto článku sme použili model strojového učenia (Rain-Net) na vykonanie zrážok nowcasting. Použili sme architektúru U-Net, ktorú sme vytvorili pomocou TensorFlow. Nahrali sme vopred pripravený model na predpovedanie satelitných snímok.

Táto implementácia sa dá vylepšiť mnohými spôsobmi. Napríklad:

Dolaďte model na svojom súbore údajov
V architektúre použite modul založený na pozornosti, ako napríklad CBAM (Konvolučný blokový modul pozornosti).

Referencie

– Blog Googlu o modeli predpovede počasia založenom na U-Net

Implementácia siete Rain-Net

– Implementácia U-Net

Príďte na jeden z našich bezplatných workshopov!

Začnite svoju kariéru ako dátový vedec s našimi bezplatnými workshopmi, ktoré sú založené na adaptabilnom učebnom pláne a vedené odborníkmi z odvetvia.