<jupyter_start><jupyter_text>Xarray Oceanography Example El Niño / Southern Oscillation (ENSO) from Sea Surface Temperature (SST)**Author:** [<NAME>](http://rabernat.github.io)According to [NOAA](https://www.esrl.noaa.gov/psd/enso/):>  El Niño and La Niña, together called the El Niño Southern Oscillation (ENSO), are periodic departures from expected sea surface temperatures (SSTs) in the equatorial Pacific Ocean. These warmer or cooler than normal ocean temperatures can affect weather patterns around the world by influencing high and low pressure systems, winds, and precipitation. ENSO may bring much needed moisture to a region while causing extremes of too much or too little water in others.In this notebook, we use the python [xarray](http://xarray.pydata.org/en/latest/) package to examine SST data from NOAA's [Extended Reconstructed Sea Surface Temperature (ERSST) v5](https://www.ncdc.noaa.gov/data-access/marineocean-data/extended-reconstructed-sea-surface-temperature-ersst-v5) product.We use [holoviews](http://holoviews.org/) to interactively visualize the data.We then demonstrate how easy it is to calculate the [Niño3.4 index](https://www.ncdc.noaa.gov/teleconnections/enso/indicators/sst/).Finally, we use the [EOFS](https://ajdawson.github.io/eofs/index.html) package to compute the Emprical Orthogoal Functions of global SST.<jupyter_code>import xarray as xr
from matplotlib import pyplot as plt
import numpy as np
import hvplot.xarray
import holoviews as hv
hv.extension('bokeh')
%matplotlib inline<jupyter_output><empty_output><jupyter_text>Access data directly from NOAA using OpenDAP protocal<jupyter_code>url = 'http://www.esrl.noaa.gov/psd/thredds/dodsC/Datasets/noaa.ersst.v5/sst.mnmean.nc'
ds = xr.open_dataset(url)
ds<jupyter_output><empty_output><jupyter_text>Select only the more recent part of the data.<jupyter_code>ds = ds.sel(time=slice('1960', '2018'))
ds
# size of data in MB
ds.nbytes/1e6 
# actually download the data
ds.load()<jupyter_output><empty_output><jupyter_text>Interactive Plot<jupyter_code>%output holomap='scrubber' fps=3
ds.sst.hvplot('lon', 'lat', cmap='Magma').redim.range(sst=(-2, 30))<jupyter_output><empty_output><jupyter_text>Calculate Climatology and Monthly AnomalyXarray makes this particularly easy<jupyter_code>sst_clim = ds.sst.groupby('time.month').mean(dim='time')
sst_anom = ds.sst.groupby('time.month') - sst_clim<jupyter_output><empty_output><jupyter_text>Detrend signal<jupyter_code>from scipy.signal import detrend
sst_anom_detrended = xr.apply_ufunc(detrend, sst_anom.fillna(0),
                                    kwargs={'axis': 0}).where(~sst_anom.isnull())<jupyter_output><empty_output><jupyter_text>Plot global mean<jupyter_code># For a global average, we need to weigh the points by cosine of latitude.
# This is not built into xarray because xarray is not specific to geoscientific data.
weights = np.cos(np.deg2rad(ds.lat)).where(~sst_anom[0].isnull())
weights /= weights.mean()

(sst_anom * weights).mean(dim=['lon', 'lat']).plot(label='raw')
(sst_anom_detrended * weights).mean(dim=['lon', 'lat']).plot(label='detrended')
plt.grid()
plt.legend()
sst_anom_detrended.hvplot('lon', 'lat', cmap='RdBu_r').redim.range(sst=(-2, 2))<jupyter_output><empty_output><jupyter_text>Calculate Oceanic Niño Index<jupyter_code>sst_anom_nino34 = sst_anom_detrended.sel(lat=slice(5, -5), lon=slice(190, 240))
sst_anom_nino34_mean = sst_anom_nino34.mean(dim=('lon', 'lat'))
oni = sst_anom_nino34_mean.rolling(time=3, center=True).mean()<jupyter_output><empty_output><jupyter_text>Plot it:<jupyter_code>oni.plot()
plt.grid()
plt.ylabel('Anomaly (dec. C)');<jupyter_output><empty_output><jupyter_text>Compare to the [official version](https://www.ncdc.noaa.gov/teleconnections/enso/indicators/sst/) from NOAA: ![NOAA Nino 3.4](https://camo.githubusercontent.com/f4fb57e583e664d9b5f1861d4b062ad0a32aeac1/68747470733a2f2f7777772e6e6364632e6e6f61612e676f762f6d6f6e69746f72696e672d636f6e74656e742f74656c65636f6e6e656374696f6e732f656c6e2d662d70672e676966)  Composite the global SST on the Niño3.4 index<jupyter_code>positive_oni = ((oni>0.5).astype('b').rolling(time=5, center=True).sum()==5)
negative_oni = ((oni<0.5).astype('b').rolling(time=5, center=True).sum()==5)
positive_oni.astype('i').plot(marker='.', label='positive')
(-negative_oni.astype('i')).plot(marker='.', label='negative')
plt.legend()
sst_anom.where(positive_oni).mean(dim='time').plot()
plt.title('SST Anomaly - Positive Niño3.4');
sst_anom.where(negative_oni).mean(dim='time').plot()
plt.title('SST Anomaly - Negative Niño3.4');<jupyter_output><empty_output><jupyter_text>Calculate Empirical Orthogonal Functions<jupyter_code>from eofs.xarray import Eof
solver = Eof(sst_anom_detrended, weights=np.sqrt(weights))
eof1 = solver.eofsAsCorrelation(neofs=1)
pc1 = solver.pcs(npcs=1, pcscaling=1)
eof1.sel(mode=0).plot()
pc1.sel(mode=0).plot()
plt.grid()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.13779767962548342, 'fraction_numerical': 0.04131894972521881, 'mean_word_length': 7.217391304347826, 'pattern_counts': {'":': 0, '<': 47, '<?xml version=': 0, '>': 48, 'https://': 6, 'lorem ipsum': 0, 'www.': 5, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '11', 'dir': 'jupyter-structured-clean-dedup', 'id': '31941', 'license': 'CC0-1.0', 'n_tokens_mistral': 1983, 'n_tokens_neox': 1840, 'n_words': 344, 'path': 'xarray_oceanography_enso.ipynb', 'repo_name': 'pangeo-data/open-source-geoscience', 'size': '4922'}