diff --git a/xarray/core/dataarray.py b/xarray/core/dataarray.py index 1b135a350d1..64f21b0eb01 100644 --- a/xarray/core/dataarray.py +++ b/xarray/core/dataarray.py @@ -2971,6 +2971,43 @@ def quantile( See Also -------- numpy.nanpercentile, pandas.Series.quantile, Dataset.quantile + + Examples + -------- + + >>> da = xr.DataArray( + ... data=[[0.7, 4.2, 9.4, 1.5], [6.5, 7.3, 2.6, 1.9]], + ... coords={"x": [7, 9], "y": [1, 1.5, 2, 2.5]}, + ... dims=("x", "y"), + ... ) + + Single quantile + >>> da.quantile(0) # or da.quantile(0, dim=...) + + array(0.7) + Coordinates: + quantile float64 0.0 + >>> da.quantile(0, dim="x") + + array([0.7, 4.2, 2.6, 1.5]) + Coordinates: + * y (y) float64 1.0 1.5 2.0 2.5 + quantile float64 0.0 + + Multiple quantiles + >>> da.quantile([0, 0.5, 1]) + + array([0.7, 3.4, 9.4]) + Coordinates: + * quantile (quantile) float64 0.0 0.5 1.0 + >>> da.quantile([0, 0.5, 1], dim="x") + + array([[0.7 , 4.2 , 2.6 , 1.5 ], + [3.6 , 5.75, 6. , 1.7 ], + [6.5 , 7.3 , 9.4 , 1.9 ]]) + Coordinates: + * y (y) float64 1.0 1.5 2.0 2.5 + * quantile (quantile) float64 0.0 0.5 1.0 """ ds = self._to_temp_dataset().quantile( diff --git a/xarray/core/dataset.py b/xarray/core/dataset.py index 1793dd2d94d..61dde6a393b 100644 --- a/xarray/core/dataset.py +++ b/xarray/core/dataset.py @@ -5116,6 +5116,48 @@ def quantile( See Also -------- numpy.nanpercentile, pandas.Series.quantile, DataArray.quantile + + Examples + -------- + + >>> ds = xr.Dataset( + ... {"a": (("x", "y"), [[0.7, 4.2, 9.4, 1.5], [6.5, 7.3, 2.6, 1.9]])}, + ... coords={"x": [7, 9], "y": [1, 1.5, 2, 2.5]}, + ... ) + + Single quantile + >>> ds.quantile(0) # or ds.quantile(0, dim=...) + + Dimensions: () + Coordinates: + quantile float64 0.0 + Data variables: + a float64 0.7 + >>> ds.quantile(0, dim="x") + + Dimensions: (y: 4) + Coordinates: + * y (y) float64 1.0 1.5 2.0 2.5 + quantile float64 0.0 + Data variables: + a (y) float64 0.7 4.2 2.6 1.5 + + Multiple quantiles + >>> ds.quantile([0, 0.5, 1]) + + Dimensions: (quantile: 3) + Coordinates: + * quantile (quantile) float64 0.0 0.5 1.0 + Data variables: + a (quantile) float64 0.7 3.4 9.4 + >>> ds.quantile([0, 0.5, 1], dim="x") + + Dimensions: (quantile: 3, y: 4) + Coordinates: + * y (y) float64 1.0 1.5 2.0 2.5 + * quantile (quantile) float64 0.0 0.5 1.0 + Data variables: + a (quantile, y) float64 0.7 4.2 2.6 1.5 3.6 ... 1.7 6.5 7.3 9.4 1.9 """ if isinstance(dim, str): diff --git a/xarray/core/groupby.py b/xarray/core/groupby.py index 7e872c74d72..cb8f6538820 100644 --- a/xarray/core/groupby.py +++ b/xarray/core/groupby.py @@ -597,6 +597,59 @@ def quantile(self, q, dim=None, interpolation="linear", keep_attrs=None): -------- numpy.nanpercentile, pandas.Series.quantile, Dataset.quantile, DataArray.quantile + + Examples + -------- + + >>> da = xr.DataArray( + ... [[1.3, 8.4, 0.7, 6.9], [0.7, 4.2, 9.4, 1.5], [6.5, 7.3, 2.6, 1.9]], + ... coords={"x": [0, 0, 1], "y": [1, 1, 2, 2]}, + ... dims=("y", "y"), + ... ) + >>> ds = xr.Dataset({"a": da}) + + Single quantile + >>> da.groupby("x").quantile(0) + + array([[0.7, 4.2, 0.7, 1.5], + [6.5, 7.3, 2.6, 1.9]]) + Coordinates: + quantile float64 0.0 + * y (y) int64 1 1 2 2 + * x (x) int64 0 1 + >>> ds.groupby("y").quantile(0, dim=...) + + Dimensions: (y: 2) + Coordinates: + quantile float64 0.0 + * y (y) int64 1 2 + Data variables: + a (y) float64 0.7 0.7 + + Multiple quantiles + >>> da.groupby("x").quantile([0, 0.5, 1]) + + array([[[0.7 , 1. , 1.3 ], + [4.2 , 6.3 , 8.4 ], + [0.7 , 5.05, 9.4 ], + [1.5 , 4.2 , 6.9 ]], + + [[6.5 , 6.5 , 6.5 ], + [7.3 , 7.3 , 7.3 ], + [2.6 , 2.6 , 2.6 ], + [1.9 , 1.9 , 1.9 ]]]) + Coordinates: + * y (y) int64 1 1 2 2 + * quantile (quantile) float64 0.0 0.5 1.0 + * x (x) int64 0 1 + >>> ds.groupby("y").quantile([0, 0.5, 1], dim=...) + + Dimensions: (quantile: 3, y: 2) + Coordinates: + * quantile (quantile) float64 0.0 0.5 1.0 + * y (y) int64 1 2 + Data variables: + a (y, quantile) float64 0.7 5.35 8.4 0.7 2.25 9.4 """ if dim is None: dim = self._group_dim