def reduce(
        self,
        func: Callable,
        dim: Dims = None,
        *,
        keep_attrs: bool | None = None,
        keepdims: bool = False,
        numeric_only: bool = False,
        **kwargs: Any,
    ) -> Self:
        """Reduce this dataset by applying `func` along some dimension(s).

        Parameters
        ----------
        func : callable
            Function which can be called in the form
            `f(x, axis=axis, **kwargs)` to return the result of reducing 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 `func` is
            applied over all dimensions.
        keep_attrs : bool or None, optional
            If True, the dataset'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.
        keepdims : bool, default: False
            If True, the dimensions which are reduced are left in the result
            as dimensions of size one. Coordinates that use these dimensions
            are removed.
        numeric_only : bool, default: False
            If True, only apply ``func`` to variables with a numeric dtype.
        **kwargs : Any
            Additional keyword arguments passed on to ``func``.

        Returns
        -------
        reduced : Dataset
            Dataset with this object's DataArrays replaced with new DataArrays
            of summarized data and the indicated dimension(s) removed.

        Examples
        --------

        >>> dataset = xr.Dataset(
        ...     {
        ...         "math_scores": (
        ...             ["student", "test"],
        ...             [[90, 85, 92], [78, 80, 85], [95, 92, 98]],
        ...         ),
        ...         "english_scores": (
        ...             ["student", "test"],
        ...             [[88, 90, 92], [75, 82, 79], [93, 96, 91]],
        ...         ),
        ...     },
        ...     coords={
        ...         "student": ["Alice", "Bob", "Charlie"],
        ...         "test": ["Test 1", "Test 2", "Test 3"],
        ...     },
        ... )

        # Calculate the 75th percentile of math scores for each student using np.percentile

        >>> percentile_scores = dataset.reduce(np.percentile, q=75, dim="test")
        >>> percentile_scores
        <xarray.Dataset> Size: 132B
        Dimensions:         (student: 3)
        Coordinates:
          * student         (student) <U7 84B 'Alice' 'Bob' 'Charlie'
        Data variables:
            math_scores     (student) float64 24B 91.0 82.5 96.5
            english_scores  (student) float64 24B 91.0 80.5 94.5
        """
        if kwargs.get("axis") is not None:
            raise ValueError(
                "passing 'axis' to Dataset reduce methods is ambiguous."
                " Please use 'dim' instead."
            )

        dims = parse_dims_as_set(dim, set(self._dims.keys()))

        if keep_attrs is None:
            keep_attrs = _get_keep_attrs(default=False)

        variables: dict[Hashable, Variable] = {}
        for name, var in self._variables.items():
            reduce_dims = [d for d in var.dims if d in dims]
            if name in self.coords:
                if not reduce_dims:
                    variables[name] = var
            elif (
                # Some reduction functions (e.g. std, var) need to run on variables
                # that don't have the reduce dims: PR5393
                not is_extension_array_dtype(var.dtype)
                and (
                    not reduce_dims
                    or not numeric_only
                    or np.issubdtype(var.dtype, np.number)
                    or (var.dtype == np.bool_)
                )
            ):
                # prefer to aggregate over axis=None rather than
                # axis=(0, 1) if they will be equivalent, because
                # the former is often more efficient
                # keep single-element dims as list, to support Hashables
                reduce_maybe_single = (
                    None
                    if len(reduce_dims) == var.ndim and var.ndim != 1
                    else reduce_dims
                )
                variables[name] = var.reduce(
                    func,
                    dim=reduce_maybe_single,
                    keep_attrs=keep_attrs,
                    keepdims=keepdims,
                    **kwargs,
                )

        coord_names = {k for k in self.coords if k in variables}
        indexes = {k: v for k, v in self._indexes.items() if k in variables}
        attrs = self.attrs if keep_attrs else None
        return self._replace_with_new_dims(
            variables, coord_names=coord_names, attrs=attrs, indexes=indexes
        )
Reduce this dataset by applying `func` along some dimension(s).

        Parameters
        ----------
        func : callable
            Function which can be called in the form
            `f(x, axis=axis, **kwargs)` to return the result of reducing 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 `func` is
            applied over all dimensions.
        keep_attrs : bool or None, optional
            If True, the dataset'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.
        keepdims : bool, default: False
            If True, the dimensions which are reduced are left in the result
            as dimensions of size one. Coordinates that use these dimensions
            are removed.
        numeric_only : bool, default: False
            If True, only apply ``func`` to variables with a numeric dtype.
        **kwargs : Any
            Additional keyword arguments passed on to ``func``.

        Returns
        -------
        reduced : Dataset
            Dataset with this object's DataArrays replaced with new DataArrays
            of summarized data and the indicated dimension(s) removed.

        Examples
        --------

        >>> dataset = xr.Dataset(
        ...     {
        ...         "math_scores": (
        ...             ["student", "test"],
        ...             [[90, 85, 92], [78, 80, 85], [95, 92, 98]],
        ...         ),
        ...         "english_scores": (
        ...             ["student", "test"],
        ...             [[88, 90, 92], [75, 82, 79], [93, 96, 91]],
        ...         ),
        ...     },
        ...     coords={
        ...         "student": ["Alice", "Bob", "Charlie"],
        ...         "test": ["Test 1", "Test 2", "Test 3"],
        ...     },
        ... )

        # Calculate the 75th percentile of math scores for each student using np.percentile

        >>> percentile_scores = dataset.reduce(np.percentile, q=75, dim="test")
        >>> percentile_scores
        <xarray.Dataset> Size: 132B
        Dimensions:         (student: 3)
        Coordinates:
          * student         (student) <U7 84B 'Alice' 'Bob' 'Charlie'
        Data variables:
            math_scores     (student) float64 24B 91.0 82.5 96.5
            english_scores  (student) float64 24B 91.0 80.5 94.5
        
xarray/core/dataset.py