Multidimensional data processing

Multidimensional data processing#

Setelah membaca chapter ini, pembaca diharapkan dapat

Memahami struktur data multidimensi (NetCDF) dan cara membacanya menggunakan Xarray.

Mengolah data multidimensi dengan operasi manipulasi dan komputasi dasar.

Intro#

Pada mata pelatihan sebelumnya, kita telah mencoba mengolah data tabular menggunakan pandas. Namun, data meteorologi-oseanografi, khususnya yang berasal dari keluaran model, tidak mudah disimpan dalam bentuk tabular.

Data keluaran model terdiri dari banyak dimensi (waktu, latitude, longitude, kedalaman/ketinggian). Misalnya, parameter temperatur pada suatu tempat di muka bumi direpresentasikan dengan tiga dimensi.

\[ T(x,y,t) \]

dimana \(x, y,\) dan \(t\) masing-masing menandakan longitude, latitude, dan waktu.

Pada bagian ini, kita akan mempelajari bagaimana xarray dapat dengan mudah untuk mengolah struktur data seperti ini.

Struktur data multi dimensi#

Setidaknya terdapat empat terminologi penting dalam data multi dimensi (yang disimpan sepeti dalam format netcdf), yaitu variabel, koordinat, dimensi atau indeks, serta atribut.

variabel: array multi dimensi berisi variabel fisis seperti suhu udara, presipitasi, tinggi gelombang, dll. Direpresentasikan oleh koordinat (latitude, longitude, dan waktu).
koordinat: wadah seperti dictionary yang menyimpan label dan array dua dimensi berisi dimensi.
dimensi: dimensi dari masing-masing koordinat (longitude, latitude, time).
atribut: dictionary untuk menyimpan metadata.

***

Lebih spesifik, struktur data utama dalam Xarray terbagi menjadi dua, yaitu DataArray dan Dataset.

DataArray#

Ini merupakan implementasi dari xarray untuk menyimpan array multi dimensi. Beberapa properti penting antara lain:

values: array multi dimensi yang menyimpan nilai dari satu variabel tertentu
dims: dimensi dari masing-masing sumbu
coords: wadah seperti dictionary yang menyimpan label dan dimensi
attrs: dictionary penyimpan metadata

Dataset#

Xarray Dataset berisi kumpulan dari DataArray yang disimpan dalam wadah berbentuk dictionary. Properti penting antara lain:

data_vars: wadah seperti dictionary untuk menyimpan DataArray
dims: dimensi dari masing-masing sumbu
coords: wadah seperti dictionary yang menyimpan label dan dimensi
attrs: dictionary penyimpan metadata

Membaca dataset#

Built-in method xarray untuk membuka dataset adalah xarray.open_dataset. Method ini dapat digunakan untuk membuka beberapa macam format data multidimensi seperti netcdf, hdf5, zarr, grib, pydap, dengan menspesifikkan argument engine.

Xarray juga mendukung untuk membaca kumpulan file dalam satu kali proses dengan method xarray.open_mfdataset. Method ini akan mereturn xarray dataset yang telah mengalami pembagian (chunking) secara otomatis.

Dalam praktik ini, kita akan mencoba untuk membuka dataset hindcast dari BMKG-OFS yang tersimpan dalam database opendap.

Load module#

# load module

import xarray as xr
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

Load data#

Multi file#

lifile = [
    "http://web.pusmar.id:8082/opendap/inawaves_hindcast/2025/H_hires_202501.nc",
    "http://web.pusmar.id:8082/opendap/inawaves_hindcast/2025/H_hires_202502.nc"
]
dsetmf = xr.open_mfdataset(lifile)
dsetmf

<xarray.Dataset> Size: 27GB
Dimensions:  (time: 474, lat: 481, lon: 881)
Coordinates:
  * time     (time) datetime64[ns] 4kB 2025-01-01 ... 2025-03-01
  * lat      (lat) float32 2kB -15.0 -14.94 -14.88 -14.81 ... 14.88 14.94 15.0
  * lon      (lon) float32 4kB 90.0 90.06 90.12 90.19 ... 144.9 144.9 145.0
Data variables: (12/17)
    hs       (time, lat, lon) float64 2GB dask.array<chunksize=(249, 481, 881), meta=np.ndarray>
    hmax     (time, lat, lon) float64 2GB dask.array<chunksize=(249, 481, 881), meta=np.ndarray>
    dir      (time, lat, lon) float64 2GB dask.array<chunksize=(249, 481, 881), meta=np.ndarray>
    dp       (time, lat, lon) float64 2GB dask.array<chunksize=(249, 481, 881), meta=np.ndarray>
    lm       (time, lat, lon) float64 2GB dask.array<chunksize=(249, 481, 881), meta=np.ndarray>
    t01      (time, lat, lon) float64 2GB dask.array<chunksize=(249, 481, 881), meta=np.ndarray>
    ...       ...
    ptp00    (time, lat, lon) float64 2GB dask.array<chunksize=(249, 481, 881), meta=np.ndarray>
    ptp01    (time, lat, lon) float64 2GB dask.array<chunksize=(249, 481, 881), meta=np.ndarray>
    ptp02    (time, lat, lon) float64 2GB dask.array<chunksize=(249, 481, 881), meta=np.ndarray>
    pdi00    (time, lat, lon) float64 2GB dask.array<chunksize=(249, 481, 881), meta=np.ndarray>
    pdi01    (time, lat, lon) float64 2GB dask.array<chunksize=(249, 481, 881), meta=np.ndarray>
    pdi02    (time, lat, lon) float64 2GB dask.array<chunksize=(249, 481, 881), meta=np.ndarray>
Attributes:
    source:       Inawaves - BMKG Ocean Forecast System (OFS)
    description:  Inawaves Model - Hindcast
    institution:  BMKG - Center For Marine Meteorology
    email:        produksi.maritim@bmkg.go.id
    Conventions:  CF-1.8

# QUICK VIEW TO THE DATA
dset['hs'][0].plot()

<matplotlib.collections.QuadMesh at 0x7fbaf82c42f0>

../_images/9aefaf0ec06f81e024f2b7ecdf209f68df7c02b852fe6c0750a4135eb19e7aac.png

Manipulasi dataset#

Operasi manipulasi dasar#

Operasi umum yang digunakan untuk mengolah dataset xarray antara lain indexing, slicing, interpolating, dan masking.

Operasi dasar		Method
Indexing		regular numpy indexing atau `.sel` atau `.isel`
Slicing		regular numpy slicing atau `.sel` atau `.isel`
Interpolating		`.sel` atau `.interp`
Masking		`.where`

Method .sel dan .isel memiliki kemampuan yang sama untuk melakukan seleksi data baik untuk DataArray maupun Dataset (regular numpy indexing hanya berlaku untuk DataArray), namun keduanya terdapat perbedaan mendasar.

	`.sel`	`.isel`
Seleksi data	berdasarkan nilai dari label dimensi	berdasarkan index atau urutan dari nilai dimensi

Apakah sadar jika terdapat langkah yang terlewat, yang ada di pemrosesan data tabular tapi tidak ada di sini? Ya, kita tidak mempertimbangkan proses pembersihan data karena output model biasanya sudah melalui proses quality control, berbeda dengan raw data tabular dari observasi. Sehingga, fokus kita dititikberatkan pada pemrosesan data lebih lanjut.

Interpolating#

Interpolasi berguna apabila kita perlu menarik data titik yang lokasinya tidak berada persis di titik grid. Selain itu, method .interp juga berguna untuk proses regriding.

Interpolasi linear diselesaikan menggunakan metode interpolasi bilinear.

***

%%time
# .sel nearest
dssel = dset['hs'].isel(time=0).sel(lat=-7.21, lon=121.21, method='nearest')
print(dssel.data)

0.791015625
CPU times: user 6.32 ms, sys: 82 μs, total: 6.4 ms
Wall time: 141 ms

%%time
# .interp
dssel = dset['hs'].isel(time=0).interp(lat=-7.21, lon=121.21, method='linear')
print(dssel.data)

0.7919578124999956
CPU times: user 145 ms, sys: 20 ms, total: 165 ms
Wall time: 2.06 s

# .interp for regridding
dssel = dset['hs'].isel(time=0).sel(lat=slice(-7.5,-6.25), lon=slice(110,118))

new_lon = np.linspace(dssel.lon[0].data, dssel.lon[-1].data, dssel.sizes['lon']*3)
new_lat = np.linspace(dssel.lat[0].data, dssel.lat[-1].data, dssel.sizes['lat']*3)
dsinterp= dssel.interp(lat=new_lat, lon=new_lon)

# Check datanya
fig, ax = plt.subplots(figsize=(20,5), ncols=2)
dssel.plot(ax=ax[0])
dsinterp.plot(ax=ax[1])
ax[0].set_title("Raw data")
ax[1].set_title("Interpolated data")
plt.tight_layout()
plt.show()

../_images/58954ceea7243473077275ef0ac36b7caaf591fe04d11de010d8976c414d370e.png

Masking#

Sebagai percobaan, kita akan memasking daratan khususnya untuk variabel hs.

# Membuat dataset unmasked

dsetfilled = dset[['hs']].isel(time=slice(0,1)).sel(lon=slice(100,120), lat=slice(-9, -6)).fillna(0.0)
dsetfilled['hs'].plot(figsize=(20,5))

<matplotlib.collections.QuadMesh at 0x7fbae97f7390>

../_images/6090453d896407c92d85124e5188312aaac8f5b18347ffc160043b783e13e23e.png

dsetnan = xr.where(dsetfilled != 0, dsetfilled['hs'], np.nan)
dsetnan['hs'].plot(figsize=(20,5))

<matplotlib.collections.QuadMesh at 0x7fbae96e8b90>

../_images/eb42405f7a587b5b68a5e6301b601135d87b960a45a852517a6197490d692529.png

Operasi komputasi dasar#

Penggabungan#

Terdapat dua method umum untuk menggabungkan dua dataset, yaitu .merge dan .concat.

`.merge`	`.concat`
Menggabungkan beberapa Dataset yang memiliki variabel yang berbeda, namun dimensi yang sama	Menggabungkan beberapa dataset yang memiliki variabel yang sama, namun dimensi yang berbeda

concat#

Bagaimana jika kita menggabungkan dataset yang memiliki variabel dan dimensi yang berbeda? Menggunakan .concat, namun variabel akan terbroadcast dengan nilai NaN apabila di salah satu dataset tidak terdapat variabel tersebut. Perhatikan contoh berikut.

Dataset yang diatas menghasilkan nilai NaN untuk variabel uwave dan vwave pada dimensi waktu kedua.

Kita bisa memeriksanya dengan method .isnull dan .sum untuk melihat jumlah dan proporsi nilainya.

df_nan_pct = (
    (dsetconcat.isnull().sum(dim=('lat', 'lon')) / (dsetconcat.sizes['lat'] * dsetconcat.sizes['lon']) * 100)
    .to_dataframe()
)
df_nan_pct

	hs	dir	uwave	vwave
time
2025-01-01 00:00:00	18.336751	18.336751	18.336751	18.336751
2025-01-01 03:00:00	18.336751	18.336751	100.000000	100.000000

dfcek = (
    (dsetconcat.isnull().sum(dim=('lat','lon')) / 
    (dsetconcat.sizes['lat'] * dsetconcat.sizes['lon']))*100
).to_dataframe()
dfcek

	hs	dir	uwave	vwave
time
2025-01-01 00:00:00	18.336751	18.336751	18.336751	18.336751
2025-01-01 03:00:00	18.336751	18.336751	100.000000	100.000000

Agregasi#

Metode agregasi yang umum digunakan untuk mereduksi dimensi Dataset antaran lain parameter statistik (.mean, .max, .min, .std, .sum, .count, .quantile, .median), .groupby, .resample, dan .rolling.

Method	Keterangan
Parameter statistik (`.mean`, `.max`, `.min`, `.std`, `.sum`, `.count`, `.quantile`, `.median`)	Agregasi berdasarkan perhitungan statistik
`.groupby`	Mengelompokkan data berdasarkan label waktu atau kategori lain
`.resample`	Mengubah frekuensi data waktu (upsample atau downsample)
`.rolling`	Menghitung agregasi dengan jendela geser (rolling window)

Agregasi statistik#

%%time

dsmean = dsetsel.mean(dim=['lon', 'lat'])
dsmax  = dsetsel.max(dim=['lon', 'lat'])
dsmin  = dsetsel.min(dim=['lon', 'lat'])
dsstd  = dsetsel.std(dim=['lon', 'lat'])
dssum  = dsetsel.sum(dim=['lon', 'lat'])
dscount  = dsetsel.count(dim=['lon', 'lat'])
dsquantile  = dsetsel.quantile(.75, dim=['lon', 'lat'])
dsmedian  = dsetsel.median(dim=['lon', 'lat'])

lids = [dsmean, dsmax, dsmin, dsstd, dssum, dscount, dsquantile, dsmedian]

CPU times: user 187 ms, sys: 8.75 ms, total: 196 ms
Wall time: 2.42 s

ncols, nrows = 4, 2
fig, ax = plt.subplots(figsize=(15, 10), ncols=ncols, nrows=nrows)
for i in range(ncols * nrows):
    col = i % ncols
    row = i // ncols
    lids[i]['hs'].plot(ax=ax[row, col])  # Note: ax[row, col] for 2D subplot
    ax[row, col].set_title(f'Plot {i}: {lids[i].attrs.get("title", "")}')

../_images/bb481e7ffe7eb6733f0ecdc85e8111839656ba852e2bc49f28d6fa78fc53bce2.png

rolling#

%%time

dsetrol = dsetsel.rolling(time=8, center=True)
dsetrol

CPU times: user 4.77 ms, sys: 25 μs, total: 4.8 ms
Wall time: 3.9 ms

DatasetRolling [time->8(center)]

# PERBANDINGAN
fig, ax = plt.subplots(figsize=(12,5), ncols=2)
dsetrol.mean()['hs'].mean(['lon','lat']).plot(ax=ax[0])
dsetsel['hs'].mean(['lon','lat']).plot(ax=ax[1])
ax[0].set_title("Rolling window (8 window mean)")
ax[1].set_title("No Rolling window")
plt.tight_layout()
plt.show()

../_images/100da98d34ebef588d826ee0c91e51d8c3f10f1fedbea407edab94d7fb144380.png

Multidimensional data processing

Contents