question about datetime_unit of xarray.DataArray.differentiate

[jesieleo]

I have a piece of data that looks like this

<xarray.Dataset>
Dimensions:    (time: 612, LEV: 15, latitude: 20, longitude: 357)
Coordinates:
  * time       (time) datetime64[ns] 1960-01-15 1960-02-15 ... 2010-12-15
  * LEV        (LEV) float64 5.01 15.07 25.28 35.76 ... 149.0 171.4 197.8 229.5
  * latitude   (latitude) float64 -4.75 -4.25 -3.75 -3.25 ... 3.75 4.25 4.75
  * longitude  (longitude) float64 114.2 114.8 115.2 115.8 ... 291.2 291.8 292.2
Data variables:
    u          (time, LEV, latitude, longitude) float32 ...
Attributes: (12/30)
    cdm_data_type:              Grid
    Conventions:                COARDS, CF-1.6, ACDD-1.3
    creator_email:              chepurin@umd.edu
    creator_name:               APDRC
    creator_type:               institution
    creator_url:                https://www.atmos.umd.edu/~ocean/
    ...                         ...
    standard_name_vocabulary:   CF Standard Name Table v29
    summary:                    Simple Ocean Data Assimilation (SODA) soda po...
    time_coverage_end:          2010-12-15T00:00:00Z
    time_coverage_start:        1983-01-15T00:00:00Z
    title:                      SODA soda pop2.2.4 [TIME][LEV][LAT][LON]
    Westernmost_Easting:        118.25

when i try to use xarray.DataArray.differentiate
data.u.differentiate('time',datetime_unit='M')
will appear

Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "D:\Anaconda3\lib\site-packages\xarray\core\dataarray.py", line 3609, in differentiate
    ds = self._to_temp_dataset().differentiate(coord, edge_order, datetime_unit)
  File "D:\Anaconda3\lib\site-packages\xarray\core\dataset.py", line 6372, in differentiate
    coord_var = coord_var._to_numeric(datetime_unit=datetime_unit)
  File "D:\Anaconda3\lib\site-packages\xarray\core\variable.py", line 2428, in _to_numeric
    numeric_array = duck_array_ops.datetime_to_numeric(
  File "D:\Anaconda3\lib\site-packages\xarray\core\duck_array_ops.py", line 466, in datetime_to_numeric
    array = array / np.timedelta64(1, datetime_unit)
TypeError: Cannot get a common metadata divisor for Numpy datatime metadata [ns] and [M] because they have incompatible nonlinear base time units.

Would you please told me is this a BUG?

[spencerkclark]

I think this is actually expected behavior despite what the documentation suggests. A datetime_unit of months ("M") or years ("Y") cannot be supported, since those values do not represent a consistent amount of time (i.e. months and years vary in their number of days). Any other values for datetime_unit should work, however, e.g. "s" or "D" etc. The underlying frequency of the data is not important; the datetime_unit only controls the units with which the time deltas are expressed in the finite difference approximation of the derivative.

We should update the documentation and potentially provide a nicer error message.

[jesieleo]

So you mean that ("M" and "Y") in datetime_unit is not applicable here, so when do these two apply? Is there any simple example?I tried to convert datetime64[ns] to datetime64[M], but failed. #5750 Since xarray converts all times to datetime64[ns], does this cause the ('Y' and 'M') functions to actually be unusable?

[spencerkclark]

That's correct. The documentation is unfortunately misleading (#8006; thank you for bringing this to our attention!) and these options do not apply, in part for the reason you stated.

You could work around this by manually converting your time coordinate to have units of ordinal months since some reference date, though it will have a somewhat tenuous physical meaning. Consider the following example:

>>> import pandas as pd; import xarray as xr
>>> times = pd.date_range("2000-05", periods=10, freq="7MS")
>>> da = xr.DataArray(range(10), dims=["time"], coords=[times])
>>> reference_date = da.time.min()
>>> months = da.time.dt.month
>>> years = da.time.dt.year
>>> ordinal_months = months + 12 * (years - reference_date.dt.year) - reference_date.dt.month
>>> da = da.assign_coords(ordinal_months=ordinal_months)
>>> da.differentiate("ordinal_months")
<xarray.DataArray (time: 10)>
array([0.14285714, 0.14285714, 0.14285714, 0.14285714, 0.14285714,
       0.14285714, 0.14285714, 0.14285714, 0.14285714, 0.14285714])
Coordinates:
  * time            (time) datetime64[ns] 2000-05-01 2000-12-01 ... 2005-08-01
    ordinal_months  (time) int64 0 7 14 21 28 35 42 49 56 63