Intraseasonal Weather Types
Contents
4. Intraseasonal Weather Types#
Obtain Intraseasonal Weather Types (IWT, at daily scale) following the methodology explained in Anderson et al. (2019)
inputs required:
Daily values of Madden-Julian Oscillation (MJO) parameters (rmm1, rmm2, phase, mjo)
in this notebook:
Obtain MJO categories (25) based on rmm1, rmm2, and phase
Fit the autoregressive logistic model with a markov order 3 and seasonality
Time-series simulation of n simulations of 1000 years of the 25 categories
Randomly obtain pairs of rmm1, rmm2 and phase from the simulated time-series
Workflow:
Intra-seasonal Weather Types (IWTs) are representative of the Madden-Julian Oscillation (MJO), which is a broad region of anomalous atmospheric circulation and convective precipitation anomalies that propagates eastward around the equator on one to two-month timescales (Madden & Julian, 1972) and exhibits correlations with relevant coastal climatology such as tropical cyclone genesis (Slade & Malony, 2013) and surface wind wave anomalies (Marshal et al., 2015, Godoi et al., 2019)
#!/usr/bin/env python
# -*- coding: utf-8 -*-
# common
import os
import os.path as op
from datetime import date, timedelta, datetime
# pip
import numpy as np
import xarray as xr
# DEV: override installed teslakit
import sys
sys.path.insert(0, op.join(os.path.abspath(''), '..','..', '..', '..'))
# teslakit
from bluemath.teslakit2.wts.mjo import categories
from bluemath.teslakit2.toolkit.alr import ALR_WRP
from bluemath.teslakit2.util.operations import GetRepeatedValues
from bluemath.teslakit2.io.aux import save_nc
from bluemath.teslakit2.plotting.mjo import Plot_MJO_phases, Plot_MJO_Categories
Warning: ecCodes 2.21.0 or higher is recommended. You are running version 2.20.0
4.1. Files and paths#
# project path
p_data = r'/media/administrador/HD2/SamoaTonga/data'
site = 'Samoa'
p_site = op.join(p_data, site)
# deliverable folder
p_deliv = op.join(p_site, 'd09_TESLA')
# output path
p_out = op.join(p_deliv,'MJO')
if not os.path.isdir(p_out): os.makedirs(p_out)
# input data
mjo_hist_file = op.join(p_deliv,'resources','MJO_hist.nc')
# output data
mjo_sim_file = op.join(p_out,'MJO_sim.nc') # simulated MJO
alr_path = op.join(p_out,'alr_w') # path to store ALR results
4.2. Parameters#
# --------------------------------------
# load data and set parameters
MJO_hist = xr.open_dataset(mjo_hist_file) # historical MJO
# MJO ALR parameters
alr_markov_order = 3
alr_seasonality = [2, 4, 8]
# Simulation
num_sims = 100
d1_sim = np.datetime64('2020-01-01').astype(datetime)
d2_sim = np.datetime64('2100-01-01').astype(datetime)
4.3. MJO phases and categories#
It is common practice in the MJO literature to separate the longitudinal location of the center of convection into eight longitudinal phases (Wheeler & Hendon, 2004). This convention was preserved in a daily index in this study, intended to be a proxy for intra-seasonal MJO oscillations by clustering the two leading PCs (henceforth referred to as IPC1 and IPC2) of outgoing longwave radiation into eight pre-defined longitudinal phases and further separating into three categories of low, medium, and high convection strength (analogous to conventions in Lafleur et al., (2015)) (Figure 4). A separate cluster was created for times when the location of the MJO is considered to have low certainty (when the vector magnitude of PC1 and PC2 is less than 1 (Wheeler & Hendon, 2004)). Altogether, the 25 clusters of Intra-seasonal Weather Types (IWTs) effectively create categorical MJO.
# --------------------------------------
# Calculate MJO categories (25 used)
rmm1 = MJO_hist['rmm1']
rmm2 = MJO_hist['rmm2']
phase = MJO_hist['phase']
categ, d_rmm_categ = categories(rmm1, rmm2, phase)
MJO_hist['categ'] = (('time',), categ)
print(MJO_hist)
<xarray.Dataset>
Dimensions: (time: 15199)
Coordinates:
* time (time) datetime64[ns] 1979-01-01 1979-01-02 ... 2020-08-11
Data variables:
mjo (time) float64 ...
phase (time) int64 6 7 7 7 7 7 7 6 6 6 7 7 7 ... 3 3 4 4 4 4 4 4 5 5 5 5
rmm1 (time) float64 0.1425 -0.2042 -0.1586 ... 1.251 0.7866 0.3321
rmm2 (time) float64 1.05 1.374 1.539 1.46 ... 0.1439 0.1006 0.2651
categ (time) int64 22 23 15 23 23 23 23 22 14 ... 12 12 12 12 13 21 25 25
# plot MJO phases
Plot_MJO_phases(rmm1, rmm2, phase);
# plot MJO categories
Plot_MJO_Categories(rmm1, rmm2, categ);
4.4. Autoregressive Logistic Regression#
Synthetic time series of the MJO are obtained with a Markov chain of the predefined IWT categorical states (statistically significant to the third order) and subsequent sampling from joint distributions of IPC1 and IPC2 within each cluster. When consecutive days in the synthetic record are sampled from the same categorical state, the randomly picked EOF pairs are ordered to preserve counterclockwise propagation of the MJO around the globe in a consistent direction.
# --------------------------------------
# Autoregressive Logistic Regression - fit model
# MJO historical data for fitting
bmus_fit = xr.Dataset(
{
'bmus' :(('time',), MJO_hist.categ.values[:]),
},
{'time' : MJO_hist.time.values[:]}
)
# ALR terms
d_terms_settings = {
'mk_order' : alr_markov_order,
'constant' : True,
'seasonality': (True, alr_seasonality),
}
# ALR wrapper
ALRW = ALR_WRP(alr_path)
ALRW.SetFitData(25, bmus_fit, d_terms_settings)
# ALR model fitting
ALRW.FitModel(max_iter=10000)
Fitting autoregressive logistic model ...
Optimization done in 30.96 seconds
# show fit report
ALRW.Report_Fit()
# --------------------------------------
# Autoregressive Logistic Regression - simulate
# simulation dates
dates_sim = [d1_sim + timedelta(days=i) for i in range((d2_sim-d1_sim).days+1)]
# launch simulation
ALR_sim = ALRW.Simulate(num_sims, dates_sim)
ALR model fit : 1979-01-01 --- 2020-08-11
ALR model sim : 2020-01-01 --- 2100-01-01
Launching 100 simulations...
# --------------------------------------
# MJO rmm1, rmm2, phase generation
# solve each ALR simulation
l_MJO_sim = []
for s in ALR_sim.n_sim:
evbmus_sim = ALR_sim.sel(n_sim=s).evbmus_sims.values[:]
# Generate rmm1 and rmm2 simulated values
rmm12_sim = np.empty((len(evbmus_sim), 2)) * np.nan
mjo_sim = np.empty(len(evbmus_sim)) * np.nan
phase_sim = np.empty(len(evbmus_sim)) * np.nan
categs = np.unique(evbmus_sim)
for c in categs:
c_ix = np.where(evbmus_sim==c)[0]
# select random values for rmm1, rmm2
options = d_rmm_categ['cat_{0}'.format(int(c))]
r = np.random.randint(options.shape[0], size=len(c_ix))
rmm12_sim[c_ix,:] = options[r,:]
# calculate mjo and phase
mjo_sim = np.sqrt(rmm12_sim[:,0]**2 + rmm12_sim[:,1]**2)
phase_sim = np.arctan2(rmm12_sim[:,0], rmm12_sim[:,1])
# internally reorder days with same category (counter-clockwise phase ascend)
l_ad = GetRepeatedValues(evbmus_sim)
for s,e in l_ad:
# get sort index by MJO phase value
ixs = np.argsort(phase_sim[s:e])
# sort mjo
rmm12_sim[s:e,0] = rmm12_sim[s:e,0][ixs]
rmm12_sim[s:e,1] = rmm12_sim[s:e,1][ixs]
mjo_sim[s:e] = mjo_sim[s:e][ixs]
phase_sim[s:e] = phase_sim[s:e][ixs]
# append simulated PCs
l_MJO_sim.append(
xr.Dataset(
{
'mjo' :(('time',), mjo_sim),
'phase' :(('time',), phase_sim),
'rmm1' :(('time',), rmm12_sim[:,0]),
'rmm2' :(('time',), rmm12_sim[:,1]),
'evbmus_sims' :(('time',), evbmus_sim),
},
{'time' : dates_sim}
)
)
# concatenate simulations
MJO_sim = xr.concat(l_MJO_sim, 'n_sim')
# store simulated MJO
save_nc(MJO_sim, mjo_sim_file, safe_time=True)
4.5. Validation#
Synthetic and historical MJO categories comparison:
Perpetual Year
Cluster Transition Probabilities
# show simulation report
ALRW.Report_Sim();
