You can do this in a single step with xarray.apply_ufunc(), which is a sort of more flexible/powerful interface to xarray's broadcasting arithmetic. Extending the toy weather example from the docs: ```python import xarray as xr import numpy as np import pandas as pd import seaborn as sns # pandas aware plotting library

np.random.seed(123)

times = pd.date_range('2000-01-01', '2001-12-31', name='time') annual_cycle = np.sin(2 * np.pi * (np.array(times.dayofyear) / 365.25 - 0.28))

base = 10 + 15 * annual_cycle.reshape(-1, 1) tmin_values = base + 3 * np.random.randn(annual_cycle.size, 3) tmax_values = base + 10 + 3 * np.random.randn(annual_cycle.size, 3)

ds = xr.Dataset({'tmin': (('time', 'location'), tmin_values), 'tmax': (('time', 'location'), tmax_values)},((62, 3), (3,), (3,)) {'time': times, 'location': ['IA', 'IN', 'IL']})

new code

ds_mean = ds.groupby('time.month').mean('time') ds_std = ds.groupby('time.month').std('time')

xarray.apply_ufunc(lambda x, m, s: (x - m) / s, ds.groupby('time.month'), ds_mean, ds_std) ```

The other way (about twice as slow) is to chain two calls to groupby(): python (ds.groupby('time.month') - ds_mean).groupby('time.month') / ds_std

I'll mark this as a documentation issue in case anyone wants to add an example to the docs.