Repository: pydata/xarray

Issue:
Different behavior of groupby with / without flox
### What happened?

(Following up from https://github.com/pydata/xarray/issues/8263)

Without flox:

```

[ins] In [1]: da =  xr.tutorial.open_dataset("air_temperature")['air']
         ...:
         ...: da.drop_vars('lat').groupby('lat').sum()
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
Cell In[1], line 3
      1 da =  xr.tutorial.open_dataset("air_temperature")['air']
----> 3 da.drop_vars('lat').groupby('lat').sum()

File ~/workspace/xarray/xarray/core/_aggregations.py:6176, in DataArrayGroupByAggregations.sum(self, dim, skipna, min_count, keep_attrs, **kwargs)
   6166     return self._flox_reduce(
   6167         func="sum",
   6168         dim=dim,
   (...)
   6173         **kwargs,
   6174     )
   6175 else:
-> 6176     return self._reduce_without_squeeze_warn(
   6177         duck_array_ops.sum,
   6178         dim=dim,
   6179         skipna=skipna,
   6180         min_count=min_count,
   6181         keep_attrs=keep_attrs,
   6182         **kwargs,
   6183     )

File ~/workspace/xarray/xarray/core/groupby.py:1417, in DataArrayGroupByBase._reduce_without_squeeze_warn(self, func, dim, axis, keep_attrs, keepdims, shortcut, **kwargs)
   1415 with warnings.catch_warnings():
   1416     warnings.filterwarnings("ignore", message="The `squeeze` kwarg")
-> 1417     check_reduce_dims(dim, self.dims)
   1419 return self._map_maybe_warn(reduce_array, shortcut=shortcut, warn_squeeze=False)

File ~/workspace/xarray/xarray/core/groupby.py:66, in check_reduce_dims(reduce_dims, dimensions)
     64     reduce_dims = [reduce_dims]
     65 if any(dim not in dimensions for dim in reduce_dims):
---> 66     raise ValueError(
     67         f"cannot reduce over dimensions {reduce_dims!r}. expected either '...' "
     68         f"to reduce over all dimensions or one or more of {dimensions!r}."
     69         f" Try passing .groupby(..., squeeze=False)"
     70     )

ValueError: cannot reduce over dimensions ['lat']. expected either '...' to reduce over all dimensions or one or more of ('time', 'lon'). Try passing .groupby(..., squeeze=False)

```

But with flox:

```
[ins] In [1]: da =  xr.tutorial.open_dataset("air_temperature")['air']
         ...:
         ...: da.drop_vars('lat').groupby('lat').sum()
Out[1]:
<xarray.DataArray 'air' (time: 2920, lat: 25, lon: 53)> Size: 31MB
array([[[241.2 , 242.5 , 243.5 , ..., 232.8 , 235.5 , 238.6 ],
        [243.8 , 244.5 , 244.7 , ..., 232.8 , 235.3 , 239.3 ],
        [250.  , 249.8 , 248.89, ..., 233.2 , 236.39, 241.7 ],
        ...,
```

### What did you expect to happen?

Identical behavior

### Minimal Complete Verifiable Example

```Python
As above
```


### MVCE confirmation

- [X] Minimal example — the example is as focused as reasonably possible to demonstrate the underlying issue in xarray.
- [X] Complete example — the example is self-contained, including all data and the text of any traceback.
- [X] Verifiable example — the example copy & pastes into an IPython prompt or [Binder notebook](https://mybinder.org/v2/gh/pydata/xarray/main?urlpath=lab/tree/doc/examples/blank_template.ipynb), returning the result.
- [X] New issue — a search of GitHub Issues suggests this is not a duplicate.
- [X] Recent environment — the issue occurs with the latest version of xarray and its dependencies.

### Relevant log output

_No response_

### Anything else we need to know?

_No response_

### Environment

<details>


INSTALLED VERSIONS
------------------
commit: d0048ef8e5c07b5290e9d2c4c144d23a45ab91f7
python: 3.11.9 (main, Apr  2 2024, 08:25:04) [Clang 15.0.0 (clang-1500.3.9.4)]
python-bits: 64
OS: Darwin
OS-release: 23.5.0
machine: arm64
processor: arm
byteorder: little
LC_ALL: en_US.UTF-8
LANG: None
LOCALE: ('en_US', 'UTF-8')
libhdf5: 1.12.2
libnetcdf: 4.9.3-development

xarray: 2024.6.1.dev50+g7b08a948
pandas: 2.2.2
numpy: 1.26.4
scipy: 1.13.1
netCDF4: 1.6.5
pydap: None
h5netcdf: 1.3.0
h5py: 3.11.0
zarr: 2.18.2
cftime: 1.6.3
nc_time_axis: 1.4.1
iris: None
bottleneck: 1.3.8
dask: 2024.7.0
distributed: 2024.7.0
matplotlib: 3.9.0
cartopy: None
seaborn: 0.13.2
numbagg: 0.8.1
fsspec: 2024.5.0
cupy: None
pint: None
sparse: None
flox: 0.9.8
numpy_groupies: 0.11.1
setuptools: 69.2.0
pip: 24.0
conda: None
pytest: 8.2.1
mypy: 1.8.0
IPython: 8.25.0
sphinx: None

</details>


Discussion:

Patch:
diff --git a/doc/whats-new.rst b/doc/whats-new.rst
index b47af127..52d587ac 100644
--- a/doc/whats-new.rst
+++ b/doc/whats-new.rst
@@ -49,10 +49,10 @@ Breaking changes
   ``origin`` or ``offset`` parameters is recommended as a replacement for using
   the ``base`` parameter and using time offset arithmetic is recommended as a
   replacement for using the ``loffset`` parameter.
-
-
-Deprecations
-~~~~~~~~~~~~
+- The ``squeeze`` kwarg to ``groupby`` is completely deprecated. This has been the source of some quite confusing
+  behaviour and has been deprecated since v2024.01.0. `groupby`` behavior is now always consistent
+  with the existing ``.groupby(..., squeeze=False)`` behavior.
+  By `Deepak Cherian <https://github.com/dcherian>`_. (:pull:`9280`)
 
 
 Bug fixes
diff --git a/xarray/core/_aggregations.py b/xarray/core/_aggregations.py
index 96f860b3..acc534d8 100644
--- a/xarray/core/_aggregations.py
+++ b/xarray/core/_aggregations.py
@@ -2316,19 +2316,6 @@ class DataArrayAggregations:
 class DatasetGroupByAggregations:
     _obj: Dataset
 
-    def _reduce_without_squeeze_warn(
-        self,
-        func: Callable[..., Any],
-        dim: Dims = None,
-        *,
-        axis: int | Sequence[int] | None = None,
-        keep_attrs: bool | None = None,
-        keepdims: bool = False,
-        shortcut: bool = True,
-        **kwargs: Any,
-    ) -> Dataset:
-        raise NotImplementedError()
-
     def reduce(
         self,
         func: Callable[..., Any],
@@ -2436,7 +2423,7 @@ class DatasetGroupByAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.count,
                 dim=dim,
                 numeric_only=False,
@@ -2532,7 +2519,7 @@ class DatasetGroupByAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.array_all,
                 dim=dim,
                 numeric_only=False,
@@ -2628,7 +2615,7 @@ class DatasetGroupByAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.array_any,
                 dim=dim,
                 numeric_only=False,
@@ -2741,7 +2728,7 @@ class DatasetGroupByAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.max,
                 dim=dim,
                 skipna=skipna,
@@ -2855,7 +2842,7 @@ class DatasetGroupByAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.min,
                 dim=dim,
                 skipna=skipna,
@@ -2971,7 +2958,7 @@ class DatasetGroupByAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.mean,
                 dim=dim,
                 skipna=skipna,
@@ -3105,7 +3092,7 @@ class DatasetGroupByAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.prod,
                 dim=dim,
                 skipna=skipna,
@@ -3240,7 +3227,7 @@ class DatasetGroupByAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.sum,
                 dim=dim,
                 skipna=skipna,
@@ -3372,7 +3359,7 @@ class DatasetGroupByAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.std,
                 dim=dim,
                 skipna=skipna,
@@ -3504,7 +3491,7 @@ class DatasetGroupByAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.var,
                 dim=dim,
                 skipna=skipna,
@@ -3606,7 +3593,7 @@ class DatasetGroupByAggregations:
         Data variables:
             da       (labels) float64 24B nan 2.0 1.5
         """
-        return self._reduce_without_squeeze_warn(
+        return self.reduce(
             duck_array_ops.median,
             dim=dim,
             skipna=skipna,
@@ -3705,7 +3692,7 @@ class DatasetGroupByAggregations:
         Data variables:
             da       (time) float64 48B 1.0 2.0 3.0 3.0 4.0 nan
         """
-        return self._reduce_without_squeeze_warn(
+        return self.reduce(
             duck_array_ops.cumsum,
             dim=dim,
             skipna=skipna,
@@ -3804,7 +3791,7 @@ class DatasetGroupByAggregations:
         Data variables:
             da       (time) float64 48B 1.0 2.0 3.0 0.0 4.0 nan
         """
-        return self._reduce_without_squeeze_warn(
+        return self.reduce(
             duck_array_ops.cumprod,
             dim=dim,
             skipna=skipna,
@@ -3817,19 +3804,6 @@ class DatasetGroupByAggregations:
 class DatasetResampleAggregations:
     _obj: Dataset
 
-    def _reduce_without_squeeze_warn(
-        self,
-        func: Callable[..., Any],
-        dim: Dims = None,
-        *,
-        axis: int | Sequence[int] | None = None,
-        keep_attrs: bool | None = None,
-        keepdims: bool = False,
-        shortcut: bool = True,
-        **kwargs: Any,
-    ) -> Dataset:
-        raise NotImplementedError()
-
     def reduce(
         self,
         func: Callable[..., Any],
@@ -3937,7 +3911,7 @@ class DatasetResampleAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.count,
                 dim=dim,
                 numeric_only=False,
@@ -4033,7 +4007,7 @@ class DatasetResampleAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.array_all,
                 dim=dim,
                 numeric_only=False,
@@ -4129,7 +4103,7 @@ class DatasetResampleAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.array_any,
                 dim=dim,
                 numeric_only=False,
@@ -4242,7 +4216,7 @@ class DatasetResampleAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.max,
                 dim=dim,
                 skipna=skipna,
@@ -4356,7 +4330,7 @@ class DatasetResampleAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.min,
                 dim=dim,
                 skipna=skipna,
@@ -4472,7 +4446,7 @@ class DatasetResampleAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.mean,
                 dim=dim,
                 skipna=skipna,
@@ -4606,7 +4580,7 @@ class DatasetResampleAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.prod,
                 dim=dim,
                 skipna=skipna,
@@ -4741,7 +4715,7 @@ class DatasetResampleAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.sum,
                 dim=dim,
                 skipna=skipna,
@@ -4873,7 +4847,7 @@ class DatasetResampleAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.std,
                 dim=dim,
                 skipna=skipna,
@@ -5005,7 +4979,7 @@ class DatasetResampleAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.var,
                 dim=dim,
                 skipna=skipna,
@@ -5107,7 +5081,7 @@ class DatasetResampleAggregations:
         Data variables:
             da       (time) float64 24B 1.0 2.0 nan
         """
-        return self._reduce_without_squeeze_warn(
+        return self.reduce(
             duck_array_ops.median,
             dim=dim,
             skipna=skipna,
@@ -5206,7 +5180,7 @@ class DatasetResampleAggregations:
         Data variables:
             da       (time) float64 48B 1.0 2.0 5.0 5.0 2.0 nan
         """
-        return self._reduce_without_squeeze_warn(
+        return self.reduce(
             duck_array_ops.cumsum,
             dim=dim,
             skipna=skipna,
@@ -5305,7 +5279,7 @@ class DatasetResampleAggregations:
         Data variables:
             da       (time) float64 48B 1.0 2.0 6.0 0.0 2.0 nan
         """
-        return self._reduce_without_squeeze_warn(
+        return self.reduce(
             duck_array_ops.cumprod,
             dim=dim,
             skipna=skipna,
@@ -5318,19 +5292,6 @@ class DatasetResampleAggregations:
 class DataArrayGroupByAggregations:
     _obj: DataArray
 
-    def _reduce_without_squeeze_warn(
-        self,
-        func: Callable[..., Any],
-        dim: Dims = None,
-        *,
-        axis: int | Sequence[int] | None = None,
-        keep_attrs: bool | None = None,
-        keepdims: bool = False,
-        shortcut: bool = True,
-        **kwargs: Any,
-    ) -> DataArray:
-        raise NotImplementedError()
-
     def reduce(
         self,
         func: Callable[..., Any],
@@ -5432,7 +5393,7 @@ class DataArrayGroupByAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.count,
                 dim=dim,
                 keep_attrs=keep_attrs,
@@ -5521,7 +5482,7 @@ class DataArrayGroupByAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.array_all,
                 dim=dim,
                 keep_attrs=keep_attrs,
@@ -5610,7 +5571,7 @@ class DataArrayGroupByAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.array_any,
                 dim=dim,
                 keep_attrs=keep_attrs,
@@ -5714,7 +5675,7 @@ class DataArrayGroupByAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.max,
                 dim=dim,
                 skipna=skipna,
@@ -5819,7 +5780,7 @@ class DataArrayGroupByAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.min,
                 dim=dim,
                 skipna=skipna,
@@ -5926,7 +5887,7 @@ class DataArrayGroupByAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.mean,
                 dim=dim,
                 skipna=skipna,
@@ -6049,7 +6010,7 @@ class DataArrayGroupByAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.prod,
                 dim=dim,
                 skipna=skipna,
@@ -6173,7 +6134,7 @@ class DataArrayGroupByAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.sum,
                 dim=dim,
                 skipna=skipna,
@@ -6294,7 +6255,7 @@ class DataArrayGroupByAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.std,
                 dim=dim,
                 skipna=skipna,
@@ -6415,7 +6376,7 @@ class DataArrayGroupByAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.var,
                 dim=dim,
                 skipna=skipna,
@@ -6509,7 +6470,7 @@ class DataArrayGroupByAggregations:
         Coordinates:
           * labels   (labels) object 24B 'a' 'b' 'c'
         """
-        return self._reduce_without_squeeze_warn(
+        return self.reduce(
             duck_array_ops.median,
             dim=dim,
             skipna=skipna,
@@ -6604,7 +6565,7 @@ class DataArrayGroupByAggregations:
           * time     (time) datetime64[ns] 48B 2001-01-31 2001-02-28 ... 2001-06-30
             labels   (time) <U1 24B 'a' 'b' 'c' 'c' 'b' 'a'
         """
-        return self._reduce_without_squeeze_warn(
+        return self.reduce(
             duck_array_ops.cumsum,
             dim=dim,
             skipna=skipna,
@@ -6699,7 +6660,7 @@ class DataArrayGroupByAggregations:
           * time     (time) datetime64[ns] 48B 2001-01-31 2001-02-28 ... 2001-06-30
             labels   (time) <U1 24B 'a' 'b' 'c' 'c' 'b' 'a'
         """
-        return self._reduce_without_squeeze_warn(
+        return self.reduce(
             duck_array_ops.cumprod,
             dim=dim,
             skipna=skipna,
@@ -6711,19 +6672,6 @@ class DataArrayGroupByAggregations:
 class DataArrayResampleAggregations:
     _obj: DataArray
 
-    def _reduce_without_squeeze_warn(
-        self,
-        func: Callable[..., Any],
-        dim: Dims = None,
-        *,
-        axis: int | Sequence[int] | None = None,
-        keep_attrs: bool | None = None,
-        keepdims: bool = False,
-        shortcut: bool = True,
-        **kwargs: Any,
-    ) -> DataArray:
-        raise NotImplementedError()
-
     def reduce(
         self,
         func: Callable[..., Any],
@@ -6825,7 +6773,7 @@ class DataArrayResampleAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.count,
                 dim=dim,
                 keep_attrs=keep_attrs,
@@ -6914,7 +6862,7 @@ class DataArrayResampleAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.array_all,
                 dim=dim,
                 keep_attrs=keep_attrs,
@@ -7003,7 +6951,7 @@ class DataArrayResampleAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.array_any,
                 dim=dim,
                 keep_attrs=keep_attrs,
@@ -7107,7 +7055,7 @@ class DataArrayResampleAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.max,
                 dim=dim,
                 skipna=skipna,
@@ -7212,7 +7160,7 @@ class DataArrayResampleAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.min,
                 dim=dim,
                 skipna=skipna,
@@ -7319,7 +7267,7 @@ class DataArrayResampleAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.mean,
                 dim=dim,
                 skipna=skipna,
@@ -7442,7 +7390,7 @@ class DataArrayResampleAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.prod,
                 dim=dim,
                 skipna=skipna,
@@ -7566,7 +7514,7 @@ class DataArrayResampleAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.sum,
                 dim=dim,
                 skipna=skipna,
@@ -7687,7 +7635,7 @@ class DataArrayResampleAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.std,
                 dim=dim,
                 skipna=skipna,
@@ -7808,7 +7756,7 @@ class DataArrayResampleAggregations:
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.var,
                 dim=dim,
                 skipna=skipna,
@@ -7902,7 +7850,7 @@ class DataArrayResampleAggregations:
         Coordinates:
           * time     (time) datetime64[ns] 24B 2001-01-31 2001-04-30 2001-07-31
         """
-        return self._reduce_without_squeeze_warn(
+        return self.reduce(
             duck_array_ops.median,
             dim=dim,
             skipna=skipna,
@@ -7997,7 +7945,7 @@ class DataArrayResampleAggregations:
             labels   (time) <U1 24B 'a' 'b' 'c' 'c' 'b' 'a'
         Dimensions without coordinates: time
         """
-        return self._reduce_without_squeeze_warn(
+        return self.reduce(
             duck_array_ops.cumsum,
             dim=dim,
             skipna=skipna,
@@ -8092,7 +8040,7 @@ class DataArrayResampleAggregations:
             labels   (time) <U1 24B 'a' 'b' 'c' 'c' 'b' 'a'
         Dimensions without coordinates: time
         """
-        return self._reduce_without_squeeze_warn(
+        return self.reduce(
             duck_array_ops.cumprod,
             dim=dim,
             skipna=skipna,
diff --git a/xarray/core/dataarray.py b/xarray/core/dataarray.py
index aef28caf..dea73cf9 100644
--- a/xarray/core/dataarray.py
+++ b/xarray/core/dataarray.py
@@ -6701,7 +6701,7 @@ class DataArray(
             Hashable | DataArray | IndexVariable | Mapping[Any, Grouper] | None
         ) = None,
         *,
-        squeeze: bool | None = None,
+        squeeze: Literal[False] = False,
         restore_coord_dims: bool = False,
         **groupers: Grouper,
     ) -> DataArrayGroupBy:
@@ -6713,10 +6713,8 @@ class DataArray(
             Array whose unique values should be used to group this array. If a
             Hashable, must be the name of a coordinate contained in this dataarray. If a dictionary,
             must map an existing variable name to a :py:class:`Grouper` instance.
-        squeeze : bool, default: True
-            If "group" is a dimension of any arrays in this dataset, `squeeze`
-            controls whether the subarrays have a dimension of length 1 along
-            that dimension or if the dimension is squeezed out.
+        squeeze : False
+            This argument is deprecated.
         restore_coord_dims : bool, default: False
             If True, also restore the dimension order of multi-dimensional
             coordinates.
@@ -6813,7 +6811,6 @@ class DataArray(
         return DataArrayGroupBy(
             self,
             (rgrouper,),
-            squeeze=squeeze,
             restore_coord_dims=restore_coord_dims,
         )
 
@@ -6826,7 +6823,7 @@ class DataArray(
         labels: ArrayLike | Literal[False] | None = None,
         precision: int = 3,
         include_lowest: bool = False,
-        squeeze: bool | None = None,
+        squeeze: Literal[False] = False,
         restore_coord_dims: bool = False,
         duplicates: Literal["raise", "drop"] = "raise",
     ) -> DataArrayGroupBy:
@@ -6858,10 +6855,8 @@ class DataArray(
             The precision at which to store and display the bins labels.
         include_lowest : bool, default: False
             Whether the first interval should be left-inclusive or not.
-        squeeze : bool, default: True
-            If "group" is a dimension of any arrays in this dataset, `squeeze`
-            controls whether the subarrays have a dimension of length 1 along
-            that dimension or if the dimension is squeezed out.
+        squeeze : False
+            This argument is deprecated.
         restore_coord_dims : bool, default: False
             If True, also restore the dimension order of multi-dimensional
             coordinates.
@@ -6909,7 +6904,6 @@ class DataArray(
         return DataArrayGroupBy(
             self,
             (rgrouper,),
-            squeeze=squeeze,
             restore_coord_dims=restore_coord_dims,
         )
 
diff --git a/xarray/core/dataset.py b/xarray/core/dataset.py
index 3d3a051c..381cf8dd 100644
--- a/xarray/core/dataset.py
+++ b/xarray/core/dataset.py
@@ -10275,7 +10275,7 @@ class Dataset(
             Hashable | DataArray | IndexVariable | Mapping[Any, Grouper] | None
         ) = None,
         *,
-        squeeze: bool | None = None,
+        squeeze: Literal[False] = False,
         restore_coord_dims: bool = False,
         **groupers: Grouper,
     ) -> DatasetGroupBy:
@@ -10354,7 +10354,6 @@ class Dataset(
         return DatasetGroupBy(
             self,
             (rgrouper,),
-            squeeze=squeeze,
             restore_coord_dims=restore_coord_dims,
         )
 
@@ -10367,7 +10366,7 @@ class Dataset(
         labels: ArrayLike | None = None,
         precision: int = 3,
         include_lowest: bool = False,
-        squeeze: bool | None = None,
+        squeeze: Literal[False] = False,
         restore_coord_dims: bool = False,
         duplicates: Literal["raise", "drop"] = "raise",
     ) -> DatasetGroupBy:
@@ -10399,10 +10398,8 @@ class Dataset(
             The precision at which to store and display the bins labels.
         include_lowest : bool, default: False
             Whether the first interval should be left-inclusive or not.
-        squeeze : bool, default: True
-            If "group" is a dimension of any arrays in this dataset, `squeeze`
-            controls whether the subarrays have a dimension of length 1 along
-            that dimension or if the dimension is squeezed out.
+        squeeze : False
+            This argument is deprecated.
         restore_coord_dims : bool, default: False
             If True, also restore the dimension order of multi-dimensional
             coordinates.
@@ -10450,7 +10447,6 @@ class Dataset(
         return DatasetGroupBy(
             self,
             (rgrouper,),
-            squeeze=squeeze,
             restore_coord_dims=restore_coord_dims,
         )
 
diff --git a/xarray/core/groupby.py b/xarray/core/groupby.py
index c335211e..ca1006b9 100644
--- a/xarray/core/groupby.py
+++ b/xarray/core/groupby.py
@@ -38,7 +38,6 @@ from xarray.core.utils import (
     FrozenMappingWarningOnValuesAccess,
     contains_only_chunked_or_numpy,
     either_dict_or_kwargs,
-    emit_user_level_warning,
     hashable,
     is_scalar,
     maybe_wrap_array,
@@ -66,31 +65,9 @@ def check_reduce_dims(reduce_dims, dimensions):
             raise ValueError(
                 f"cannot reduce over dimensions {reduce_dims!r}. expected either '...' "
                 f"to reduce over all dimensions or one or more of {dimensions!r}."
-                f" Try passing .groupby(..., squeeze=False)"
             )
 
 
-def _maybe_squeeze_indices(
-    indices, squeeze: bool | None, grouper: ResolvedGrouper, warn: bool
-):
-    from xarray.core.groupers import UniqueGrouper
-
-    is_unique_grouper = isinstance(grouper.grouper, UniqueGrouper)
-    can_squeeze = is_unique_grouper and grouper.grouper.can_squeeze
-    if squeeze in [None, True] and can_squeeze:
-        if isinstance(indices, slice):
-            if indices.stop - indices.start == 1:
-                if (squeeze is None and warn) or squeeze is True:
-                    emit_user_level_warning(
-                        "The `squeeze` kwarg to GroupBy is being removed."
-                        "Pass .groupby(..., squeeze=False) to disable squeezing,"
-                        " which is the new default, and to silence this warning."
-                    )
-
-                indices = indices.start
-    return indices
-
-
 def _codes_to_group_indices(inverse: np.ndarray, N: int) -> GroupIndices:
     assert inverse.ndim == 1
     groups: GroupIndices = tuple([] for _ in range(N))
@@ -369,17 +346,15 @@ class ResolvedGrouper(Generic[T_DataWithCoords]):
             self.unique_coord = encoded.unique_coord
 
 
-def _validate_groupby_squeeze(squeeze: bool | None) -> None:
+def _validate_groupby_squeeze(squeeze: Literal[False]) -> None:
     # While we don't generally check the type of every arg, passing
     # multiple dimensions as multiple arguments is common enough, and the
     # consequences hidden enough (strings evaluate as true) to warrant
     # checking here.
     # A future version could make squeeze kwarg only, but would face
     # backward-compat issues.
-    if squeeze is not None and not isinstance(squeeze, bool):
-        raise TypeError(
-            f"`squeeze` must be None,  True or False, but {squeeze} was supplied"
-        )
+    if squeeze is not False:
+        raise TypeError(f"`squeeze` must be False, but {squeeze} was supplied.")
 
 
 def _resolve_group(
@@ -466,7 +441,6 @@ class GroupBy(Generic[T_Xarray]):
         "_unique_coord",
         "_dims",
         "_sizes",
-        "_squeeze",
         # Save unstacked object for flox
         "_original_obj",
         "_original_group",
@@ -475,7 +449,6 @@ class GroupBy(Generic[T_Xarray]):
     )
     _obj: T_Xarray
     groupers: tuple[ResolvedGrouper]
-    _squeeze: bool | None
     _restore_coord_dims: bool
 
     _original_obj: T_Xarray
@@ -492,7 +465,6 @@ class GroupBy(Generic[T_Xarray]):
         self,
         obj: T_Xarray,
         groupers: tuple[ResolvedGrouper],
-        squeeze: bool | None = False,
         restore_coord_dims: bool = True,
     ) -> None:
         """Create a GroupBy object
@@ -517,7 +489,6 @@ class GroupBy(Generic[T_Xarray]):
         # specification for the groupby operation
         self._obj = grouper.stacked_obj
         self._restore_coord_dims = restore_coord_dims
-        self._squeeze = squeeze
 
         # These should generalize to multiple groupers
         self._group_indices = grouper.group_indices
@@ -542,14 +513,8 @@ class GroupBy(Generic[T_Xarray]):
         """
         if self._sizes is None:
             (grouper,) = self.groupers
-            index = _maybe_squeeze_indices(
-                self._group_indices[0],
-                self._squeeze,
-                grouper,
-                warn=True,
-            )
+            index = self._group_indices[0]
             self._sizes = self._obj.isel({self._group_dim: index}).sizes
-
         return self._sizes
 
     def map(
@@ -582,11 +547,7 @@ class GroupBy(Generic[T_Xarray]):
         # provided to mimic pandas.groupby
         if self._groups is None:
             (grouper,) = self.groupers
-            squeezed_indices = (
-                _maybe_squeeze_indices(ind, self._squeeze, grouper, warn=idx > 0)
-                for idx, ind in enumerate(self._group_indices)
-            )
-            self._groups = dict(zip(grouper.unique_coord.values, squeezed_indices))
+            self._groups = dict(zip(grouper.unique_coord.values, self._group_indices))
         return self._groups
 
     def __getitem__(self, key: GroupKey) -> T_Xarray:
@@ -594,10 +555,7 @@ class GroupBy(Generic[T_Xarray]):
         Get DataArray or Dataset corresponding to a particular group label.
         """
         (grouper,) = self.groupers
-        index = _maybe_squeeze_indices(
-            self.groups[key], self._squeeze, grouper, warn=True
-        )
-        return self._obj.isel({self._group_dim: index})
+        return self._obj.isel({self._group_dim: self.groups[key]})
 
     def __len__(self) -> int:
         (grouper,) = self.groupers
@@ -616,13 +574,10 @@ class GroupBy(Generic[T_Xarray]):
             ", ".join(format_array_flat(grouper.full_index, 30).split()),
         )
 
-    def _iter_grouped(self, warn_squeeze=True) -> Iterator[T_Xarray]:
+    def _iter_grouped(self) -> Iterator[T_Xarray]:
         """Iterate over each element in this group"""
         (grouper,) = self.groupers
         for idx, indices in enumerate(self._group_indices):
-            indices = _maybe_squeeze_indices(
-                indices, self._squeeze, grouper, warn=warn_squeeze and idx == 0
-            )
             yield self._obj.isel({self._group_dim: indices})
 
     def _infer_concat_args(self, applied_example):
@@ -809,14 +764,6 @@ class GroupBy(Generic[T_Xarray]):
                 # set explicitly to avoid unnecessarily accumulating count
                 kwargs["min_count"] = 0
 
-        # weird backcompat
-        # reducing along a unique indexed dimension with squeeze=True
-        # should raise an error
-        if (dim is None or dim == name) and name in obj.xindexes:
-            index = obj.indexes[name]
-            if index.is_unique and self._squeeze:
-                raise ValueError(f"cannot reduce over dimensions {name!r}")
-
         unindexed_dims: tuple[Hashable, ...] = tuple()
         if isinstance(grouper.group, _DummyGroup) and not isbin:
             unindexed_dims = (name,)
@@ -1141,23 +1088,17 @@ class DataArrayGroupByBase(GroupBy["DataArray"], DataArrayGroupbyArithmetic):
     def dims(self) -> tuple[Hashable, ...]:
         if self._dims is None:
             (grouper,) = self.groupers
-            index = _maybe_squeeze_indices(
-                self._group_indices[0], self._squeeze, grouper, warn=True
-            )
+            index = self._group_indices[0]
             self._dims = self._obj.isel({self._group_dim: index}).dims
-
         return self._dims
 
-    def _iter_grouped_shortcut(self, warn_squeeze=True):
+    def _iter_grouped_shortcut(self):
         """Fast version of `_iter_grouped` that yields Variables without
         metadata
         """
         var = self._obj.variable
         (grouper,) = self.groupers
         for idx, indices in enumerate(self._group_indices):
-            indices = _maybe_squeeze_indices(
-                indices, self._squeeze, grouper, warn=warn_squeeze and idx == 0
-            )
             yield var[{self._group_dim: indices}]
 
     def _concat_shortcut(self, applied, dim, positions=None):
@@ -1237,24 +1178,7 @@ class DataArrayGroupByBase(GroupBy["DataArray"], DataArrayGroupbyArithmetic):
         applied : DataArray
             The result of splitting, applying and combining this array.
         """
-        return self._map_maybe_warn(
-            func, args, warn_squeeze=True, shortcut=shortcut, **kwargs
-        )
-
-    def _map_maybe_warn(
-        self,
-        func: Callable[..., DataArray],
-        args: tuple[Any, ...] = (),
-        *,
-        warn_squeeze: bool = True,
-        shortcut: bool | None = None,
-        **kwargs: Any,
-    ) -> DataArray:
-        grouped = (
-            self._iter_grouped_shortcut(warn_squeeze)
-            if shortcut
-            else self._iter_grouped(warn_squeeze)
-        )
+        grouped = self._iter_grouped_shortcut() if shortcut else self._iter_grouped()
         applied = (maybe_wrap_array(arr, func(arr, *args, **kwargs)) for arr in grouped)
         return self._combine(applied, shortcut=shortcut)
 
@@ -1356,68 +1280,6 @@ class DataArrayGroupByBase(GroupBy["DataArray"], DataArrayGroupbyArithmetic):
 
         return self.map(reduce_array, shortcut=shortcut)
 
-    def _reduce_without_squeeze_warn(
-        self,
-        func: Callable[..., Any],
-        dim: Dims = None,
-        *,
-        axis: int | Sequence[int] | None = None,
-        keep_attrs: bool | None = None,
-        keepdims: bool = False,
-        shortcut: bool = True,
-        **kwargs: Any,
-    ) -> DataArray:
-        """Reduce the items in this group by applying `func` along some
-        dimension(s).
-
-        Parameters
-        ----------
-        func : callable
-            Function which can be called in the form
-            `func(x, axis=axis, **kwargs)` to return the result of collapsing
-            an np.ndarray over an integer valued axis.
-        dim : "...", str, Iterable of Hashable or None, optional
-            Dimension(s) over which to apply `func`. If None, apply over the
-            groupby dimension, if "..." apply over all dimensions.
-        axis : int or sequence of int, optional
-            Axis(es) over which to apply `func`. Only one of the 'dimension'
-            and 'axis' arguments can be supplied. If neither are supplied, then
-            `func` is calculated over all dimension for each group item.
-        keep_attrs : bool, optional
-            If True, the datasets's attributes (`attrs`) will be copied from
-            the original object to the new one.  If False (default), the new
-            object will be returned without attributes.
-        **kwargs : dict
-            Additional keyword arguments passed on to `func`.
-
-        Returns
-        -------
-        reduced : Array
-            Array with summarized data and the indicated dimension(s)
-            removed.
-        """
-        if dim is None:
-            dim = [self._group_dim]
-
-        if keep_attrs is None:
-            keep_attrs = _get_keep_attrs(default=True)
-
-        def reduce_array(ar: DataArray) -> DataArray:
-            return ar.reduce(
-                func=func,
-                dim=dim,
-                axis=axis,
-                keep_attrs=keep_attrs,
-                keepdims=keepdims,
-                **kwargs,
-            )
-
-        with warnings.catch_warnings():
-            warnings.filterwarnings("ignore", message="The `squeeze` kwarg")
-            check_reduce_dims(dim, self.dims)
-
-        return self._map_maybe_warn(reduce_array, shortcut=shortcut, warn_squeeze=False)
-
 
 # https://github.com/python/mypy/issues/9031
 class DataArrayGroupBy(  # type: ignore[misc]
@@ -1436,12 +1298,7 @@ class DatasetGroupByBase(GroupBy["Dataset"], DatasetGroupbyArithmetic):
     def dims(self) -> Frozen[Hashable, int]:
         if self._dims is None:
             (grouper,) = self.groupers
-            index = _maybe_squeeze_indices(
-                self._group_indices[0],
-                self._squeeze,
-                grouper,
-                warn=True,
-            )
+            index = self._group_indices[0]
             self._dims = self._obj.isel({self._group_dim: index}).dims
 
         return FrozenMappingWarningOnValuesAccess(self._dims)
@@ -1482,18 +1339,8 @@ class DatasetGroupByBase(GroupBy["Dataset"], DatasetGroupbyArithmetic):
         applied : Dataset
             The result of splitting, applying and combining this dataset.
         """
-        return self._map_maybe_warn(func, args, shortcut, warn_squeeze=True, **kwargs)
-
-    def _map_maybe_warn(
-        self,
-        func: Callable[..., Dataset],
-        args: tuple[Any, ...] = (),
-        shortcut: bool | None = None,
-        warn_squeeze: bool = False,
-        **kwargs: Any,
-    ) -> Dataset:
         # ignore shortcut if set (for now)
-        applied = (func(ds, *args, **kwargs) for ds in self._iter_grouped(warn_squeeze))
+        applied = (func(ds, *args, **kwargs) for ds in self._iter_grouped())
         return self._combine(applied)
 
     def apply(self, func, args=(), shortcut=None, **kwargs):
@@ -1588,68 +1435,6 @@ class DatasetGroupByBase(GroupBy["Dataset"], DatasetGroupbyArithmetic):
 
         return self.map(reduce_dataset)
 
-    def _reduce_without_squeeze_warn(
-        self,
-        func: Callable[..., Any],
-        dim: Dims = None,
-        *,
-        axis: int | Sequence[int] | None = None,
-        keep_attrs: bool | None = None,
-        keepdims: bool = False,
-        shortcut: bool = True,
-        **kwargs: Any,
-    ) -> Dataset:
-        """Reduce the items in this group by applying `func` along some
-        dimension(s).
-
-        Parameters
-        ----------
-        func : callable
-            Function which can be called in the form
-            `func(x, axis=axis, **kwargs)` to return the result of collapsing
-            an np.ndarray over an integer valued axis.
-        dim : ..., str, Iterable of Hashable or None, optional
-            Dimension(s) over which to apply `func`. By default apply over the
-            groupby dimension, with "..." apply over all dimensions.
-        axis : int or sequence of int, optional
-            Axis(es) over which to apply `func`. Only one of the 'dimension'
-            and 'axis' arguments can be supplied. If neither are supplied, then
-            `func` is calculated over all dimension for each group item.
-        keep_attrs : bool, optional
-            If True, the datasets's attributes (`attrs`) will be copied from
-            the original object to the new one.  If False (default), the new
-            object will be returned without attributes.
-        **kwargs : dict
-            Additional keyword arguments passed on to `func`.
-
-        Returns
-        -------
-        reduced : Dataset
-            Array with summarized data and the indicated dimension(s)
-            removed.
-        """
-        if dim is None:
-            dim = [self._group_dim]
-
-        if keep_attrs is None:
-            keep_attrs = _get_keep_attrs(default=True)
-
-        def reduce_dataset(ds: Dataset) -> Dataset:
-            return ds.reduce(
-                func=func,
-                dim=dim,
-                axis=axis,
-                keep_attrs=keep_attrs,
-                keepdims=keepdims,
-                **kwargs,
-            )
-
-        with warnings.catch_warnings():
-            warnings.filterwarnings("ignore", message="The `squeeze` kwarg")
-            check_reduce_dims(dim, self.dims)
-
-        return self._map_maybe_warn(reduce_dataset, warn_squeeze=False)
-
     def assign(self, **kwargs: Any) -> Dataset:
         """Assign data variables by group.
 
diff --git a/xarray/core/groupers.py b/xarray/core/groupers.py
index 1b96988b..8f9810e8 100644
--- a/xarray/core/groupers.py
+++ b/xarray/core/groupers.py
@@ -74,17 +74,6 @@ class EncodedGroups:
 class Grouper(ABC):
     """Abstract base class for Grouper objects that allow specializing GroupBy instructions."""
 
-    @property
-    def can_squeeze(self) -> bool:
-        """
-        Do not use.
-
-        .. deprecated:: 2024.01.0
-            This is a deprecated method. It will be deleted when the `squeeze` kwarg is deprecated.
-            Only ``UniqueGrouper`` should override it.
-        """
-        return False
-
     @abstractmethod
     def factorize(self, group: T_Group) -> EncodedGroups:
         """
@@ -118,15 +107,6 @@ class UniqueGrouper(Grouper):
 
     _group_as_index: pd.Index | None = field(default=None, repr=False)
 
-    @property
-    def is_unique_and_monotonic(self) -> bool:
-        if isinstance(self.group, _DummyGroup):
-            return True
-        index = self.group_as_index
-        return index.is_unique and (
-            index.is_monotonic_increasing or index.is_monotonic_decreasing
-        )
-
     @property
     def group_as_index(self) -> pd.Index:
         """Caches the group DataArray as a pandas Index."""
@@ -134,16 +114,18 @@ class UniqueGrouper(Grouper):
             self._group_as_index = self.group.to_index()
         return self._group_as_index
 
-    @property
-    def can_squeeze(self) -> bool:
-        """This is a deprecated method and will be removed eventually."""
-        is_dimension = self.group.dims == (self.group.name,)
-        return is_dimension and self.is_unique_and_monotonic
-
     def factorize(self, group1d: T_Group) -> EncodedGroups:
         self.group = group1d
 
-        if self.can_squeeze:
+        index = self.group_as_index
+        is_unique_and_monotonic = isinstance(self.group, _DummyGroup) or (
+            index.is_unique
+            and (index.is_monotonic_increasing or index.is_monotonic_decreasing)
+        )
+        is_dimension = self.group.dims == (self.group.name,)
+        can_squeeze = is_dimension and is_unique_and_monotonic
+
+        if can_squeeze:
             return self._factorize_dummy()
         else:
             return self._factorize_unique()
diff --git a/xarray/core/resample.py b/xarray/core/resample.py
index ec86f2a2..a175b284 100644
--- a/xarray/core/resample.py
+++ b/xarray/core/resample.py
@@ -286,21 +286,9 @@ class DataArrayResample(Resample["DataArray"], DataArrayGroupByBase, DataArrayRe
         applied : DataArray
             The result of splitting, applying and combining this array.
         """
-        return self._map_maybe_warn(func, args, shortcut, warn_squeeze=True, **kwargs)
-
-    def _map_maybe_warn(
-        self,
-        func: Callable[..., Any],
-        args: tuple[Any, ...] = (),
-        shortcut: bool | None = False,
-        warn_squeeze: bool = True,
-        **kwargs: Any,
-    ) -> DataArray:
         # TODO: the argument order for Resample doesn't match that for its parent,
         # GroupBy
-        combined = super()._map_maybe_warn(
-            func, shortcut=shortcut, args=args, warn_squeeze=warn_squeeze, **kwargs
-        )
+        combined = super().map(func, shortcut=shortcut, args=args, **kwargs)
 
         # If the aggregation function didn't drop the original resampling
         # dimension, then we need to do so before we can rename the proxy
@@ -380,18 +368,8 @@ class DatasetResample(Resample["Dataset"], DatasetGroupByBase, DatasetResampleAg
         applied : Dataset
             The result of splitting, applying and combining this dataset.
         """
-        return self._map_maybe_warn(func, args, shortcut, warn_squeeze=True, **kwargs)
-
-    def _map_maybe_warn(
-        self,
-        func: Callable[..., Any],
-        args: tuple[Any, ...] = (),
-        shortcut: bool | None = None,
-        warn_squeeze: bool = True,
-        **kwargs: Any,
-    ) -> Dataset:
         # ignore shortcut if set (for now)
-        applied = (func(ds, *args, **kwargs) for ds in self._iter_grouped(warn_squeeze))
+        applied = (func(ds, *args, **kwargs) for ds in self._iter_grouped())
         combined = self._combine(applied)
 
         # If the aggregation function didn't drop the original resampling
@@ -466,27 +444,6 @@ class DatasetResample(Resample["Dataset"], DatasetGroupByBase, DatasetResampleAg
             **kwargs,
         )
 
-    def _reduce_without_squeeze_warn(
-        self,
-        func: Callable[..., Any],
-        dim: Dims = None,
-        *,
-        axis: int | Sequence[int] | None = None,
-        keep_attrs: bool | None = None,
-        keepdims: bool = False,
-        shortcut: bool = True,
-        **kwargs: Any,
-    ) -> Dataset:
-        return super()._reduce_without_squeeze_warn(
-            func=func,
-            dim=dim,
-            axis=axis,
-            keep_attrs=keep_attrs,
-            keepdims=keepdims,
-            shortcut=shortcut,
-            **kwargs,
-        )
-
     def asfreq(self) -> Dataset:
         """Return values of original object at the new up-sampling frequency;
         essentially a re-index with new times set to NaN.
diff --git a/xarray/util/generate_aggregations.py b/xarray/util/generate_aggregations.py
index b59dc36c..d93c94b1 100644
--- a/xarray/util/generate_aggregations.py
+++ b/xarray/util/generate_aggregations.py
@@ -91,19 +91,6 @@ GROUPBY_PREAMBLE = """
 class {obj}{cls}Aggregations:
     _obj: {obj}
 
-    def _reduce_without_squeeze_warn(
-        self,
-        func: Callable[..., Any],
-        dim: Dims = None,
-        *,
-        axis: int | Sequence[int] | None = None,
-        keep_attrs: bool | None = None,
-        keepdims: bool = False,
-        shortcut: bool = True,
-        **kwargs: Any,
-    ) -> {obj}:
-        raise NotImplementedError()
-
     def reduce(
         self,
         func: Callable[..., Any],
@@ -128,19 +115,6 @@ RESAMPLE_PREAMBLE = """
 class {obj}{cls}Aggregations:
     _obj: {obj}
 
-    def _reduce_without_squeeze_warn(
-        self,
-        func: Callable[..., Any],
-        dim: Dims = None,
-        *,
-        axis: int | Sequence[int] | None = None,
-        keep_attrs: bool | None = None,
-        keepdims: bool = False,
-        shortcut: bool = True,
-        **kwargs: Any,
-    ) -> {obj}:
-        raise NotImplementedError()
-
     def reduce(
         self,
         func: Callable[..., Any],
@@ -457,7 +431,7 @@ class GroupByAggregationGenerator(AggregationGenerator):
 
         if method_is_not_flox_supported:
             return f"""\
-        return self._reduce_without_squeeze_warn(
+        return self.reduce(
             duck_array_ops.{method.array_method},
             dim=dim,{extra_kwargs}
             keep_attrs=keep_attrs,
@@ -484,7 +458,7 @@ class GroupByAggregationGenerator(AggregationGenerator):
                 **kwargs,
             )
         else:
-            return self._reduce_without_squeeze_warn(
+            return self.reduce(
                 duck_array_ops.{method.array_method},
                 dim=dim,{extra_kwargs}
                 keep_attrs=keep_attrs,


Test patch:
diff --git a/xarray/tests/test_computation.py b/xarray/tests/test_computation.py
index b000de31..8b480b02 100644
--- a/xarray/tests/test_computation.py
+++ b/xarray/tests/test_computation.py
@@ -118,10 +118,8 @@ def test_apply_identity() -> None:
     assert_identical(variable, apply_identity(variable))
     assert_identical(data_array, apply_identity(data_array))
     assert_identical(data_array, apply_identity(data_array.groupby("x")))
-    assert_identical(data_array, apply_identity(data_array.groupby("x", squeeze=False)))
     assert_identical(dataset, apply_identity(dataset))
     assert_identical(dataset, apply_identity(dataset.groupby("x")))
-    assert_identical(dataset, apply_identity(dataset.groupby("x", squeeze=False)))
 
 
 def add(a, b):
@@ -521,10 +519,8 @@ def test_apply_output_core_dimension() -> None:
     assert_identical(stacked_variable, stack_negative(variable))
     assert_identical(stacked_data_array, stack_negative(data_array))
     assert_identical(stacked_dataset, stack_negative(dataset))
-    with pytest.warns(UserWarning, match="The `squeeze` kwarg"):
-        assert_identical(stacked_data_array, stack_negative(data_array.groupby("x")))
-    with pytest.warns(UserWarning, match="The `squeeze` kwarg"):
-        assert_identical(stacked_dataset, stack_negative(dataset.groupby("x")))
+    assert_identical(stacked_data_array, stack_negative(data_array.groupby("x")))
+    assert_identical(stacked_dataset, stack_negative(dataset.groupby("x")))
 
     def original_and_stack_negative(obj):
         def func(x):
@@ -551,13 +547,11 @@ def test_apply_output_core_dimension() -> None:
     assert_identical(dataset, out0)
     assert_identical(stacked_dataset, out1)
 
-    with pytest.warns(UserWarning, match="The `squeeze` kwarg"):
-        out0, out1 = original_and_stack_negative(data_array.groupby("x"))
+    out0, out1 = original_and_stack_negative(data_array.groupby("x"))
     assert_identical(data_array, out0)
     assert_identical(stacked_data_array, out1)
 
-    with pytest.warns(UserWarning, match="The `squeeze` kwarg"):
-        out0, out1 = original_and_stack_negative(dataset.groupby("x"))
+    out0, out1 = original_and_stack_negative(dataset.groupby("x"))
     assert_identical(dataset, out0)
     assert_identical(stacked_dataset, out1)
 
diff --git a/xarray/tests/test_groupby.py b/xarray/tests/test_groupby.py
index 9519a555..fc074060 100644
--- a/xarray/tests/test_groupby.py
+++ b/xarray/tests/test_groupby.py
@@ -62,60 +62,44 @@ def test_consolidate_slices() -> None:
 
 
 @pytest.mark.filterwarnings("ignore:return type")
-def test_groupby_dims_property(dataset, recwarn) -> None:
-    # dims is sensitive to squeeze, always warn
-    with pytest.warns(UserWarning, match="The `squeeze` kwarg"):
-        assert dataset.groupby("x").dims == dataset.isel(x=1).dims
-        assert dataset.groupby("y").dims == dataset.isel(y=1).dims
-    # in pytest-8, pytest.warns() no longer clears all warnings
-    recwarn.clear()
-
-    # when squeeze=False, no warning should be raised
-    assert tuple(dataset.groupby("x", squeeze=False).dims) == tuple(
-        dataset.isel(x=slice(1, 2)).dims
-    )
-    assert tuple(dataset.groupby("y", squeeze=False).dims) == tuple(
-        dataset.isel(y=slice(1, 2)).dims
-    )
-    assert len(recwarn) == 0
+def test_groupby_dims_property(dataset) -> None:
+    with pytest.warns(FutureWarning, match="The return type of"):
+        assert dataset.groupby("x").dims == dataset.isel(x=[1]).dims
+    with pytest.warns(FutureWarning, match="The return type of"):
+        assert dataset.groupby("y").dims == dataset.isel(y=[1]).dims
+
+    assert tuple(dataset.groupby("x").dims) == tuple(dataset.isel(x=slice(1, 2)).dims)
+    assert tuple(dataset.groupby("y").dims) == tuple(dataset.isel(y=slice(1, 2)).dims)
 
     dataset = dataset.drop_vars(["cat"])
     stacked = dataset.stack({"xy": ("x", "y")})
-    assert tuple(stacked.groupby("xy", squeeze=False).dims) == tuple(
-        stacked.isel(xy=[0]).dims
-    )
-    assert len(recwarn) == 0
+    assert tuple(stacked.groupby("xy").dims) == tuple(stacked.isel(xy=[0]).dims)
 
 
 def test_groupby_sizes_property(dataset) -> None:
-    with pytest.warns(UserWarning, match="The `squeeze` kwarg"):
-        assert dataset.groupby("x").sizes == dataset.isel(x=1).sizes
-    with pytest.warns(UserWarning, match="The `squeeze` kwarg"):
-        assert dataset.groupby("y").sizes == dataset.isel(y=1).sizes
+    assert dataset.groupby("x").sizes == dataset.isel(x=[1]).sizes
+    assert dataset.groupby("y").sizes == dataset.isel(y=[1]).sizes
     dataset = dataset.drop_vars("cat")
     stacked = dataset.stack({"xy": ("x", "y")})
-    with pytest.warns(UserWarning, match="The `squeeze` kwarg"):
-        assert stacked.groupby("xy").sizes == stacked.isel(xy=0).sizes
+    assert stacked.groupby("xy").sizes == stacked.isel(xy=[0]).sizes
 
 
 def test_multi_index_groupby_map(dataset) -> None:
     # regression test for GH873
     ds = dataset.isel(z=1, drop=True)[["foo"]]
     expected = 2 * ds
-    # The function in `map` may be sensitive to squeeze, always warn
-    with pytest.warns(UserWarning, match="The `squeeze` kwarg"):
-        actual = (
-            ds.stack(space=["x", "y"])
-            .groupby("space")
-            .map(lambda x: 2 * x)
-            .unstack("space")
-        )
+    actual = (
+        ds.stack(space=["x", "y"])
+        .groupby("space")
+        .map(lambda x: 2 * x)
+        .unstack("space")
+    )
     assert_equal(expected, actual)
 
 
 @pytest.mark.parametrize("grouper", [dict(group="x"), dict(x=UniqueGrouper())])
 def test_reduce_numeric_only(dataset, grouper: dict) -> None:
-    gb = dataset.groupby(**grouper, squeeze=False)
+    gb = dataset.groupby(**grouper)
     with xr.set_options(use_flox=False):
         expected = gb.sum()
     with xr.set_options(use_flox=True):
@@ -222,8 +206,7 @@ def test_da_groupby_map_func_args() -> None:
 
     array = xr.DataArray([1, 1, 1], [("x", [1, 2, 3])])
     expected = xr.DataArray([3, 3, 3], [("x", [1, 2, 3])])
-    with pytest.warns(UserWarning, match="The `squeeze` kwarg"):
-        actual = array.groupby("x").map(func, args=(1,), arg3=1)
+    actual = array.groupby("x").map(func, args=(1,), arg3=1)
     assert_identical(expected, actual)
 
 
@@ -233,9 +216,7 @@ def test_ds_groupby_map_func_args() -> None:
 
     dataset = xr.Dataset({"foo": ("x", [1, 1, 1])}, {"x": [1, 2, 3]})
     expected = xr.Dataset({"foo": ("x", [3, 3, 3])}, {"x": [1, 2, 3]})
-    # The function in `map` may be sensitive to squeeze, always warn
-    with pytest.warns(UserWarning, match="The `squeeze` kwarg"):
-        actual = dataset.groupby("x").map(func, args=(1,), arg3=1)
+    actual = dataset.groupby("x").map(func, args=(1,), arg3=1)
     assert_identical(expected, actual)
 
 
@@ -549,10 +530,8 @@ def test_da_groupby_assign_coords() -> None:
     actual = xr.DataArray(
         [[3, 4, 5], [6, 7, 8]], dims=["y", "x"], coords={"y": range(2), "x": range(3)}
     )
-    with pytest.warns(UserWarning, match="The `squeeze` kwarg"):
-        actual1 = actual.groupby("x").assign_coords({"y": [-1, -2]})
-    with pytest.warns(UserWarning, match="The `squeeze` kwarg"):
-        actual2 = actual.groupby("x").assign_coords(y=[-1, -2])
+    actual1 = actual.groupby("x").assign_coords({"y": [-1, -2]})
+    actual2 = actual.groupby("x").assign_coords(y=[-1, -2])
     expected = xr.DataArray(
         [[3, 4, 5], [6, 7, 8]], dims=["y", "x"], coords={"y": [-1, -2], "x": range(3)}
     )
@@ -708,11 +687,10 @@ def test_groupby_reduce_dimension_error(array) -> None:
     with pytest.raises(ValueError, match=r"cannot reduce over dimensions"):
         grouped.mean(("x", "y", "asd"))
 
-    with pytest.warns(UserWarning, match="The `squeeze` kwarg"):
-        assert_identical(array.mean("x"), grouped.reduce(np.mean, "x"))
-        assert_allclose(array.mean(["x", "z"]), grouped.reduce(np.mean, ["x", "z"]))
+    assert_identical(array.mean("x"), grouped.reduce(np.mean, "x"))
+    assert_allclose(array.mean(["x", "z"]), grouped.reduce(np.mean, ["x", "z"]))
 
-    grouped = array.groupby("y", squeeze=False)
+    grouped = array.groupby("y")
     assert_identical(array, grouped.mean())
 
     assert_identical(array.mean("x"), grouped.reduce(np.mean, "x"))
@@ -753,34 +731,19 @@ def test_groupby_none_group_name() -> None:
 
 
 def test_groupby_getitem(dataset) -> None:
-    with pytest.warns(UserWarning, match="The `squeeze` kwarg"):
-        assert_identical(dataset.sel(x="a"), dataset.groupby("x")["a"])
-    with pytest.warns(UserWarning, match="The `squeeze` kwarg"):
-        assert_identical(dataset.sel(z=1), dataset.groupby("z")[1])
-    with pytest.warns(UserWarning, match="The `squeeze` kwarg"):
-        assert_identical(dataset.foo.sel(x="a"), dataset.foo.groupby("x")["a"])
-    with pytest.warns(UserWarning, match="The `squeeze` kwarg"):
-        assert_identical(dataset.foo.sel(z=1), dataset.foo.groupby("z")[1])
-    with pytest.warns(UserWarning, match="The `squeeze` kwarg"):
-        assert_identical(dataset.cat.sel(y=1), dataset.cat.groupby("y")[1])
-
-    assert_identical(dataset.sel(x=["a"]), dataset.groupby("x", squeeze=False)["a"])
-    assert_identical(dataset.sel(z=[1]), dataset.groupby("z", squeeze=False)[1])
-
-    assert_identical(
-        dataset.foo.sel(x=["a"]), dataset.foo.groupby("x", squeeze=False)["a"]
-    )
-    assert_identical(dataset.foo.sel(z=[1]), dataset.foo.groupby("z", squeeze=False)[1])
 
-    assert_identical(dataset.cat.sel(y=[1]), dataset.cat.groupby("y", squeeze=False)[1])
+    assert_identical(dataset.sel(x=["a"]), dataset.groupby("x")["a"])
+    assert_identical(dataset.sel(z=[1]), dataset.groupby("z")[1])
+    assert_identical(dataset.foo.sel(x=["a"]), dataset.foo.groupby("x")["a"])
+    assert_identical(dataset.foo.sel(z=[1]), dataset.foo.groupby("z")[1])
+    assert_identical(dataset.cat.sel(y=[1]), dataset.cat.groupby("y")[1])
+
     with pytest.raises(
         NotImplementedError, match="Cannot broadcast 1d-only pandas categorical array."
     ):
-        dataset.groupby("boo", squeeze=False)
+        dataset.groupby("boo")
     dataset = dataset.drop_vars(["cat"])
-    actual = (
-        dataset.groupby("boo", squeeze=False)["f"].unstack().transpose("x", "y", "z")
-    )
+    actual = dataset.groupby("boo")["f"].unstack().transpose("x", "y", "z")
     expected = dataset.sel(y=[1], z=[1, 2]).transpose("x", "y", "z")
     assert_identical(expected, actual)
 
@@ -790,7 +753,7 @@ def test_groupby_dataset() -> None:
         {"z": (["x", "y"], np.random.randn(3, 5))},
         {"x": ("x", list("abc")), "c": ("x", [0, 1, 0]), "y": range(5)},
     )
-    groupby = data.groupby("x", squeeze=False)
+    groupby = data.groupby("x")
     assert len(groupby) == 3
     expected_groups = {"a": slice(0, 1), "b": slice(1, 2), "c": slice(2, 3)}
     assert groupby.groups == expected_groups
@@ -807,55 +770,25 @@ def test_groupby_dataset() -> None:
         return x
 
     for k in ["x", "c", "y"]:
-        actual2 = data.groupby(k, squeeze=False).map(identity)
+        actual2 = data.groupby(k).map(identity)
         assert_equal(data, actual2)
 
 
-def test_groupby_dataset_squeeze_None() -> None:
-    """Delete when removing squeeze."""
-    data = Dataset(
-        {"z": (["x", "y"], np.random.randn(3, 5))},
-        {"x": ("x", list("abc")), "c": ("x", [0, 1, 0]), "y": range(5)},
-    )
-    groupby = data.groupby("x")
-    assert len(groupby) == 3
-    expected_groups = {"a": 0, "b": 1, "c": 2}
-    with pytest.warns(UserWarning, match="The `squeeze` kwarg"):
-        assert groupby.groups == expected_groups
-    expected_items = [
-        ("a", data.isel(x=0)),
-        ("b", data.isel(x=1)),
-        ("c", data.isel(x=2)),
-    ]
-    with pytest.warns(UserWarning, match="The `squeeze` kwarg"):
-        for actual1, expected1 in zip(groupby, expected_items):
-            assert actual1[0] == expected1[0]
-            assert_equal(actual1[1], expected1[1])
-
-    def identity(x):
-        return x
-
-    with pytest.warns(UserWarning, match="The `squeeze` kwarg"):
-        for k in ["x", "c"]:
-            actual2 = data.groupby(k).map(identity)
-            assert_equal(data, actual2)
-
-
 def test_groupby_dataset_returns_new_type() -> None:
     data = Dataset({"z": (["x", "y"], np.random.randn(3, 5))})
 
-    actual1 = data.groupby("x", squeeze=False).map(lambda ds: ds["z"])
+    actual1 = data.groupby("x").map(lambda ds: ds["z"])
     expected1 = data["z"]
     assert_identical(expected1, actual1)
 
-    actual2 = data["z"].groupby("x", squeeze=False).map(lambda x: x.to_dataset())
+    actual2 = data["z"].groupby("x").map(lambda x: x.to_dataset())
     expected2 = data
     assert_identical(expected2, actual2)
 
 
 def test_groupby_dataset_iter() -> None:
     data = create_test_data()
-    for n, (t, sub) in enumerate(list(data.groupby("dim1", squeeze=False))[:3]):
+    for n, (t, sub) in enumerate(list(data.groupby("dim1"))[:3]):
         assert data["dim1"][n] == t
         assert_equal(data["var1"][[n]], sub["var1"])
         assert_equal(data["var2"][[n]], sub["var2"])
@@ -911,14 +844,13 @@ def test_groupby_dataset_reduce_ellipsis(by_func) -> None:
     assert_allclose(expected, actual)
 
 
-@pytest.mark.parametrize("squeeze", [True, False])
-def test_groupby_dataset_math(squeeze: bool) -> None:
+def test_groupby_dataset_math() -> None:
     def reorder_dims(x):
         return x.transpose("dim1", "dim2", "dim3", "time")
 
     ds = create_test_data()
     ds["dim1"] = ds["dim1"]
-    grouped = ds.groupby("dim1", squeeze=squeeze)
+    grouped = ds.groupby("dim1")
 
     expected = reorder_dims(ds + ds.coords["dim1"])
     actual = grouped + ds.coords["dim1"]
@@ -1027,7 +959,7 @@ def test_groupby_bins_cut_kwargs(use_flox: bool) -> None:
 
     with xr.set_options(use_flox=use_flox):
         actual = da.groupby_bins(
-            "x", bins=x_bins, include_lowest=True, right=False, squeeze=False
+            "x", bins=x_bins, include_lowest=True, right=False
         ).mean()
     expected = xr.DataArray(
         np.array([[1.0, 2.0], [5.0, 6.0], [9.0, 10.0]]),
@@ -1267,11 +1199,10 @@ class TestDataArrayGroupBy:
 
     def test_groupby_iter(self) -> None:
         for (act_x, act_dv), (exp_x, exp_ds) in zip(
-            self.dv.groupby("y", squeeze=False), self.ds.groupby("y", squeeze=False)
+            self.dv.groupby("y"), self.ds.groupby("y")
         ):
             assert exp_x == act_x
             assert_identical(exp_ds["foo"], act_dv)
-        with pytest.warns(UserWarning, match="The `squeeze` kwarg"):
             for (_, exp_dv), (_, act_dv) in zip(
                 self.dv.groupby("x"), self.dv.groupby("x")
             ):
@@ -1296,8 +1227,7 @@ class TestDataArrayGroupBy:
         "by, use_da", [("x", False), ("y", False), ("y", True), ("abc", False)]
     )
     @pytest.mark.parametrize("shortcut", [True, False])
-    @pytest.mark.parametrize("squeeze", [None, True, False])
-    def test_groupby_map_identity(self, by, use_da, shortcut, squeeze, recwarn) -> None:
+    def test_groupby_map_identity(self, by, use_da, shortcut) -> None:
         expected = self.da
         if use_da:
             by = expected.coords[by]
@@ -1305,14 +1235,10 @@ class TestDataArrayGroupBy:
         def identity(x):
             return x
 
-        grouped = expected.groupby(by, squeeze=squeeze)
+        grouped = expected.groupby(by)
         actual = grouped.map(identity, shortcut=shortcut)
         assert_identical(expected, actual)
 
-        # abc is not a dim coordinate so no warnings expected!
-        if (by.name if use_da else by) != "abc":
-            assert len(recwarn) == (1 if squeeze in [None, True] else 0)
-
     def test_groupby_sum(self) -> None:
         array = self.da
         grouped = array.groupby("abc")
@@ -1472,10 +1398,9 @@ class TestDataArrayGroupBy:
         expected = change_metadata(expected)
         assert_equal(expected, actual)
 
-    @pytest.mark.parametrize("squeeze", [True, False])
-    def test_groupby_math_squeeze(self, squeeze: bool) -> None:
+    def test_groupby_math_squeeze(self) -> None:
         array = self.da
-        grouped = array.groupby("x", squeeze=squeeze)
+        grouped = array.groupby("x")
 
         expected = array + array.coords["x"]
         actual = grouped + array.coords["x"]
@@ -1545,7 +1470,7 @@ class TestDataArrayGroupBy:
             ("a", ("a", "y")),
             ("b", ("x", "b")),
         ]:
-            result = array.groupby(by, squeeze=False).map(lambda x: x.squeeze())
+            result = array.groupby(by).map(lambda x: x.squeeze())
             assert result.dims == expected_dims
 
     def test_groupby_restore_coord_dims(self) -> None:
@@ -1565,7 +1490,7 @@ class TestDataArrayGroupBy:
             ("a", ("a", "y")),
             ("b", ("x", "b")),
         ]:
-            result = array.groupby(by, squeeze=False, restore_coord_dims=True).map(
+            result = array.groupby(by, restore_coord_dims=True).map(
                 lambda x: x.squeeze()
             )["c"]
             assert result.dims == expected_dims
@@ -2556,8 +2481,8 @@ def test_groupby_dim_no_dim_equal(use_flox: bool) -> None:
         data=[1, 2, 3, 4], dims="lat", coords={"lat": np.linspace(0, 1.01, 4)}
     )
     with xr.set_options(use_flox=use_flox):
-        actual1 = da.drop_vars("lat").groupby("lat", squeeze=False).sum()
-        actual2 = da.groupby("lat", squeeze=False).sum()
+        actual1 = da.drop_vars("lat").groupby("lat").sum()
+        actual2 = da.groupby("lat").sum()
     assert_identical(actual1, actual2.drop_vars("lat"))