Repository: pydata/xarray

Issue:
Numeric scalar variable attributes (including fill_value, scale_factor, add_offset) are 1-d instead of 0-d with h5netcdf engine, triggering ValueError: non-broadcastable output on application when loading single elements
<!-- Please include a self-contained copy-pastable example that generates the issue if possible.

Please be concise with code posted. See guidelines below on how to provide a good bug report:

- Craft Minimal Bug Reports: http://matthewrocklin.com/blog/work/2018/02/28/minimal-bug-reports
- Minimal Complete Verifiable Examples: https://stackoverflow.com/help/mcve

Bug reports that follow these guidelines are easier to diagnose, and so are often handled much more quickly.
-->

**What happened**:

When I try to open a NetCDF file using the `h5netcdf` engine, accessing a single data point before scale factors have been applied results in `ValueError: non-broadcastable output operand with shape () doesn't match the broadcast shape (1,)`.  The MCVE (see below) results in:

```
Traceback (most recent call last):
  File "mwe93.py", line 4, in <module>
    ds["Rad"][400, 300].load()
  File "/data/gholl/miniconda3/envs/py38b/lib/python3.8/site-packages/xarray/core/dataarray.py", line 808, in load
    ds = self._to_temp_dataset().load(**kwargs)
  File "/data/gholl/miniconda3/envs/py38b/lib/python3.8/site-packages/xarray/core/dataset.py", line 662, in load
    v.load()
  File "/data/gholl/miniconda3/envs/py38b/lib/python3.8/site-packages/xarray/core/variable.py", line 439, in load
    self._data = np.asarray(self._data)
  File "/data/gholl/miniconda3/envs/py38b/lib/python3.8/site-packages/numpy/core/_asarray.py", line 83, in asarray
    return array(a, dtype, copy=False, order=order)
  File "/data/gholl/miniconda3/envs/py38b/lib/python3.8/site-packages/xarray/core/indexing.py", line 685, in __array__
    self._ensure_cached()
  File "/data/gholl/miniconda3/envs/py38b/lib/python3.8/site-packages/xarray/core/indexing.py", line 682, in _ensure_cached
    self.array = NumpyIndexingAdapter(np.asarray(self.array))
  File "/data/gholl/miniconda3/envs/py38b/lib/python3.8/site-packages/numpy/core/_asarray.py", line 83, in asarray
    return array(a, dtype, copy=False, order=order)
  File "/data/gholl/miniconda3/envs/py38b/lib/python3.8/site-packages/xarray/core/indexing.py", line 655, in __array__
    return np.asarray(self.array, dtype=dtype)
  File "/data/gholl/miniconda3/envs/py38b/lib/python3.8/site-packages/numpy/core/_asarray.py", line 83, in asarray
    return array(a, dtype, copy=False, order=order)
  File "/data/gholl/miniconda3/envs/py38b/lib/python3.8/site-packages/xarray/core/indexing.py", line 560, in __array__
    return np.asarray(array[self.key], dtype=None)
  File "/data/gholl/miniconda3/envs/py38b/lib/python3.8/site-packages/numpy/core/_asarray.py", line 83, in asarray
    return array(a, dtype, copy=False, order=order)
  File "/data/gholl/miniconda3/envs/py38b/lib/python3.8/site-packages/xarray/coding/variables.py", line 70, in __array__
    return self.func(self.array)
  File "/data/gholl/miniconda3/envs/py38b/lib/python3.8/site-packages/xarray/coding/variables.py", line 220, in _scale_offset_decoding
    data *= scale_factor
ValueError: non-broadcastable output operand with shape () doesn't match the broadcast shape (1,)
```


**What you expected to happen**:

I expect the data access to work similarly as when opening with other engines.

**Minimal Complete Verifiable Example**:

```python
import xarray
fn = "/data/gholl/cache/fogtools/abi/2017/03/14/20/06/7/OR_ABI-L1b-RadF-M3C07_G16_s20170732006100_e20170732016478_c20170732016514.nc"
with xarray.open_dataset(fn, engine="h5netcdf") as ds:
    ds["Rad"][400, 300].load()
```

**Anything else we need to know?**:

An earlier version of this issue, and some comments, refer to fsspec or working on open files, but that proved to have nothing to do with the problem.  

**Environment**:

I've confirmed this issue installing xarray from latest master, which means xarray 0.16.2.dev11+gf821fe20 at the time of writing,

<details><summary>Output of <tt>xr.show_versions()</tt></summary>


INSTALLED VERSIONS
------------------
commit: None
python: 3.8.5 | packaged by conda-forge | (default, Sep 24 2020, 16:55:52) 
[GCC 7.5.0]
python-bits: 64
OS: Linux
OS-release: 4.12.14-lp150.12.82-default
machine: x86_64
processor: x86_64
byteorder: little
LC_ALL: None
LANG: en_GB.UTF-8
LOCALE: en_GB.UTF-8
libhdf5: 1.10.6
libnetcdf: 4.7.4

xarray: 0.16.2.dev11+gf821fe20
pandas: 1.1.2
numpy: 1.19.1
scipy: 1.5.2
netCDF4: 1.5.4
pydap: None
h5netcdf: 0.8.1
h5py: 2.10.0
Nio: None
zarr: 2.4.0
cftime: 1.2.1
nc_time_axis: None
PseudoNetCDF: None
rasterio: 1.1.6
cfgrib: None
iris: None
bottleneck: None
dask: 2.27.0
distributed: 2.27.0
matplotlib: 3.3.2
cartopy: None
seaborn: None
numbagg: None
pint: None
setuptools: 49.6.0.post20200917
pip: 20.2.3
conda: None
pytest: 6.0.2
IPython: None
sphinx: None

</details>


Discussion:
gerritholl: Is this bugfix notable enough to need a `whats-new.rst` entry?

For the unit test, I tried to construct an object that would emulate what is produced when reading a NetCDF4 file with the h5netcdf engine, but I gave up and settled for a temporary file instead.  If this is an undesired approach, I could use some guidance in how to construct the appropriate object that will expose the problem.
gerritholl: If this makes more sense as an integration test than as a unit test (for which I need help, see other comment), should I mark the current test in some way and/or move it to a different source file?

Patch:
diff --git a/ci/requirements/doc.yml b/ci/requirements/doc.yml
index 4170569f..6b7e022c 100644
--- a/ci/requirements/doc.yml
+++ b/ci/requirements/doc.yml
@@ -13,7 +13,7 @@ dependencies:
   - ipython
   - iris>=2.3
   - jupyter_client
-  - matplotlib-base=3.3.0
+  - matplotlib-base
   - nbsphinx
   - netcdf4>=1.5
   - numba
diff --git a/ci/requirements/py36.yml b/ci/requirements/py36.yml
index 0200a7a8..7c06bcb8 100644
--- a/ci/requirements/py36.yml
+++ b/ci/requirements/py36.yml
@@ -32,7 +32,7 @@ dependencies:
   - pip
   - pseudonetcdf
   - pydap
-  - pynio
+  # - pynio: not compatible with netCDF4>1.5.3; only tested in py36-bare-minimum
   - pytest
   - pytest-cov
   - pytest-env
diff --git a/ci/requirements/py37.yml b/ci/requirements/py37.yml
index 943ab755..c9702deb 100644
--- a/ci/requirements/py37.yml
+++ b/ci/requirements/py37.yml
@@ -32,7 +32,7 @@ dependencies:
   - pip
   - pseudonetcdf
   - pydap
-  - pynio
+  # - pynio: not compatible with netCDF4>1.5.3; only tested in py36-bare-minimum
   - pytest
   - pytest-cov
   - pytest-env
diff --git a/ci/requirements/py38-all-but-dask.yml b/ci/requirements/py38-all-but-dask.yml
index 843d762b..547418ad 100644
--- a/ci/requirements/py38-all-but-dask.yml
+++ b/ci/requirements/py38-all-but-dask.yml
@@ -29,7 +29,7 @@ dependencies:
   - pip
   - pseudonetcdf
   - pydap
-  - pynio
+  # - pynio: not compatible with netCDF4>1.5.3; only tested in py36-bare-minimum
   - pytest
   - pytest-cov
   - pytest-env
diff --git a/ci/requirements/py38.yml b/ci/requirements/py38.yml
index 7dff3a1b..b80b4fde 100644
--- a/ci/requirements/py38.yml
+++ b/ci/requirements/py38.yml
@@ -32,7 +32,7 @@ dependencies:
   - pip
   - pseudonetcdf
   - pydap
-  - pynio
+  # - pynio: not compatible with netCDF4>1.5.3; only tested in py36-bare-minimum
   - pytest
   - pytest-cov
   - pytest-env
diff --git a/doc/installing.rst b/doc/installing.rst
index 62e026e2..52820f8c 100644
--- a/doc/installing.rst
+++ b/doc/installing.rst
@@ -25,8 +25,9 @@ For netCDF and IO
 - `pydap <http://www.pydap.org/>`__: used as a fallback for accessing OPeNDAP
 - `h5netcdf <https://github.com/shoyer/h5netcdf>`__: an alternative library for
   reading and writing netCDF4 files that does not use the netCDF-C libraries
-- `pynio <https://www.pyngl.ucar.edu/Nio.shtml>`__: for reading GRIB and other
-  geoscience specific file formats. Note that pynio is not available for Windows.
+- `PyNIO <https://www.pyngl.ucar.edu/Nio.shtml>`__: for reading GRIB and other
+  geoscience specific file formats. Note that PyNIO is not available for Windows and
+  that the PyNIO backend may be moved outside of xarray in the future.
 - `zarr <http://zarr.readthedocs.io/>`__: for chunked, compressed, N-dimensional arrays.
 - `cftime <https://unidata.github.io/cftime>`__: recommended if you
   want to encode/decode datetimes for non-standard calendars or dates before
diff --git a/doc/io.rst b/doc/io.rst
index 956d9394..8010c726 100644
--- a/doc/io.rst
+++ b/doc/io.rst
@@ -105,6 +105,12 @@ Dataset and DataArray objects, and no array values are loaded into memory until
 you try to perform some sort of actual computation. For an example of how these
 lazy arrays work, see the OPeNDAP section below.
 
+There may be minor differences in the :py:class:`Dataset` object returned
+when reading a NetCDF file with different engines. For example,
+single-valued attributes are returned as scalars by the default
+``engine=netcdf4``, but as arrays of size ``(1,)`` when reading with
+``engine=h5netcdf``.
+
 It is important to note that when you modify values of a Dataset, even one
 linked to files on disk, only the in-memory copy you are manipulating in xarray
 is modified: the original file on disk is never touched.
@@ -1022,6 +1028,11 @@ We recommend installing PyNIO via conda::
 
     conda install -c conda-forge pynio
 
+    .. note::
+    
+    PyNIO is no longer actively maintained and conflicts with netcdf4 > 1.5.3.
+    The PyNIO backend may be moved outside of xarray in the future.
+
 .. _PyNIO: https://www.pyngl.ucar.edu/Nio.shtml
 
 .. _io.PseudoNetCDF:
diff --git a/doc/whats-new.rst b/doc/whats-new.rst
index 6b944fd1..448ab135 100644
--- a/doc/whats-new.rst
+++ b/doc/whats-new.rst
@@ -34,6 +34,8 @@ New Features
 Bug fixes
 ~~~~~~~~~
 
+- Fix bug where reading a scalar value from a NetCDF file opened with the ``h5netcdf`` backend would raise a ``ValueError`` when ``decode_cf=True`` (:issue:`4471`, :pull:`4485`).
+  By `Gerrit Holl <https://github.com/gerritholl>`_.
 - Fix bug where datetime64 times are silently changed to incorrect values if they are outside the valid date range for ns precision when provided in some other units (:issue:`4427`, :pull:`4454`).
   By `Andrew Pauling <https://github.com/andrewpauling>`_
 - Fix silently overwriting the `engine` key when passing :py:func:`open_dataset` a file object
diff --git a/xarray/backends/api.py b/xarray/backends/api.py
index 434ea5a1..0317c61a 100644
--- a/xarray/backends/api.py
+++ b/xarray/backends/api.py
@@ -25,7 +25,7 @@ from ..core.combine import (
     combine_by_coords,
 )
 from ..core.dataarray import DataArray
-from ..core.dataset import Dataset
+from ..core.dataset import Dataset, _maybe_chunk
 from ..core.utils import close_on_error, is_grib_path, is_remote_uri
 from .common import AbstractDataStore, ArrayWriter
 from .locks import _get_scheduler
@@ -191,14 +191,15 @@ def _validate_dataset_names(dataset):
         if isinstance(name, str):
             if not name:
                 raise ValueError(
-                    "Invalid name for DataArray or Dataset key: "
+                    f"Invalid name {name!r} for DataArray or Dataset key: "
                     "string must be length 1 or greater for "
                     "serialization to netCDF files"
                 )
         elif name is not None:
             raise TypeError(
-                "DataArray.name or Dataset key must be either a "
-                "string or None for serialization to netCDF files"
+                f"Invalid name {name!r} for DataArray or Dataset key: "
+                "must be either a string or None for serialization to netCDF "
+                "files"
             )
 
     for k in dataset.variables:
@@ -214,22 +215,22 @@ def _validate_attrs(dataset):
         if isinstance(name, str):
             if not name:
                 raise ValueError(
-                    "Invalid name for attr: string must be "
+                    f"Invalid name for attr {name!r}: string must be "
                     "length 1 or greater for serialization to "
                     "netCDF files"
                 )
         else:
             raise TypeError(
-                "Invalid name for attr: {} must be a string for "
-                "serialization to netCDF files".format(name)
+                f"Invalid name for attr: {name!r} must be a string for "
+                "serialization to netCDF files"
             )
 
         if not isinstance(value, (str, Number, np.ndarray, np.number, list, tuple)):
             raise TypeError(
-                "Invalid value for attr: {} must be a number, "
+                f"Invalid value for attr {name!r}: {value!r} must be a number, "
                 "a string, an ndarray or a list/tuple of "
                 "numbers/strings for serialization to netCDF "
-                "files".format(value)
+                "files"
             )
 
     # Check attrs on the dataset itself
@@ -457,9 +458,8 @@ def open_dataset(
 
     if backend_kwargs is None:
         backend_kwargs = {}
-    extra_kwargs = {}
 
-    def maybe_decode_store(store, chunks, lock=False):
+    def maybe_decode_store(store, chunks):
         ds = conventions.decode_cf(
             store,
             mask_and_scale=mask_and_scale,
@@ -525,7 +525,12 @@ def open_dataset(
                 chunks = dict.fromkeys(ds.dims, chunks)
 
             variables = {
-                k: store.maybe_chunk(k, v, chunks, overwrite_encoded_chunks)
+                k: _maybe_chunk(
+                    k,
+                    v,
+                    store.get_chunk(k, v, chunks),
+                    overwrite_encoded_chunks=overwrite_encoded_chunks,
+                )
                 for k, v in ds.variables.items()
             }
             ds2 = ds._replace(variables)
@@ -543,18 +548,17 @@ def open_dataset(
         if engine is None:
             engine = _autodetect_engine(filename_or_obj)
 
-        if engine in ["netcdf4", "h5netcdf"]:
+        extra_kwargs = {}
+        if group is not None:
             extra_kwargs["group"] = group
+        if lock is not None:
             extra_kwargs["lock"] = lock
-        elif engine in ["pynio", "pseudonetcdf", "cfgrib"]:
-            extra_kwargs["lock"] = lock
-        elif engine == "zarr":
+
+        if engine == "zarr":
             backend_kwargs = backend_kwargs.copy()
             overwrite_encoded_chunks = backend_kwargs.pop(
                 "overwrite_encoded_chunks", None
             )
-            extra_kwargs["mode"] = "r"
-            extra_kwargs["group"] = group
 
         opener = _get_backend_cls(engine)
         store = opener(filename_or_obj, **extra_kwargs, **backend_kwargs)
diff --git a/xarray/backends/zarr.py b/xarray/backends/zarr.py
index 260d27fb..8482a1ba 100644
--- a/xarray/backends/zarr.py
+++ b/xarray/backends/zarr.py
@@ -390,22 +390,6 @@ class ZarrStore(AbstractWritableDataStore):
                 chunk_spec[dim] = chunks[dim]
         return chunk_spec
 
-    def maybe_chunk(self, name, var, chunks, overwrite_encoded_chunks):
-        chunk_spec = self.get_chunk(name, var, chunks)
-
-        if (var.ndim > 0) and (chunk_spec is not None):
-            from dask.base import tokenize
-
-            # does this cause any data to be read?
-            token2 = tokenize(name, var._data, chunks)
-            name2 = f"xarray-{name}-{token2}"
-            var = var.chunk(chunk_spec, name=name2, lock=None)
-            if overwrite_encoded_chunks and var.chunks is not None:
-                var.encoding["chunks"] = tuple(x[0] for x in var.chunks)
-            return var
-        else:
-            return var
-
     def store(
         self,
         variables,
diff --git a/xarray/coding/variables.py b/xarray/coding/variables.py
index 92e45484..80277e9c 100644
--- a/xarray/coding/variables.py
+++ b/xarray/coding/variables.py
@@ -269,6 +269,10 @@ class CFScaleOffsetCoder(VariableCoder):
             scale_factor = pop_to(attrs, encoding, "scale_factor", name=name)
             add_offset = pop_to(attrs, encoding, "add_offset", name=name)
             dtype = _choose_float_dtype(data.dtype, "add_offset" in attrs)
+            if np.ndim(scale_factor) > 0:
+                scale_factor = scale_factor.item()
+            if np.ndim(add_offset) > 0:
+                add_offset = add_offset.item()
             transform = partial(
                 _scale_offset_decoding,
                 scale_factor=scale_factor,
diff --git a/xarray/core/dataset.py b/xarray/core/dataset.py
index 1777ee35..126c3dec 100644
--- a/xarray/core/dataset.py
+++ b/xarray/core/dataset.py
@@ -359,6 +359,34 @@ def _assert_empty(args: tuple, msg: str = "%s") -> None:
         raise ValueError(msg % args)
 
 
+def _maybe_chunk(
+    name,
+    var,
+    chunks=None,
+    token=None,
+    lock=None,
+    name_prefix="xarray-",
+    overwrite_encoded_chunks=False,
+):
+    from dask.base import tokenize
+
+    if chunks is not None:
+        chunks = {dim: chunks[dim] for dim in var.dims if dim in chunks}
+    if var.ndim:
+        # when rechunking by different amounts, make sure dask names change
+        # by provinding chunks as an input to tokenize.
+        # subtle bugs result otherwise. see GH3350
+        token2 = tokenize(name, token if token else var._data, chunks)
+        name2 = f"{name_prefix}{name}-{token2}"
+        var = var.chunk(chunks, name=name2, lock=lock)
+
+        if overwrite_encoded_chunks and var.chunks is not None:
+            var.encoding["chunks"] = tuple(x[0] for x in var.chunks)
+        return var
+    else:
+        return var
+
+
 def as_dataset(obj: Any) -> "Dataset":
     """Cast the given object to a Dataset.
 
@@ -1761,7 +1789,6 @@ class Dataset(Mapping, ImplementsDatasetReduce, DataWithCoords):
         -------
         chunked : xarray.Dataset
         """
-        from dask.base import tokenize
 
         if isinstance(chunks, (Number, str)):
             chunks = dict.fromkeys(self.dims, chunks)
@@ -1774,26 +1801,10 @@ class Dataset(Mapping, ImplementsDatasetReduce, DataWithCoords):
                     "object: %s" % bad_dims
                 )
 
-        def selkeys(dict_, keys):
-            if dict_ is None:
-                return None
-            return {d: dict_[d] for d in keys if d in dict_}
-
-        def maybe_chunk(name, var, chunks):
-            chunks = selkeys(chunks, var.dims)
-            if not chunks:
-                chunks = None
-            if var.ndim > 0:
-                # when rechunking by different amounts, make sure dask names change
-                # by provinding chunks as an input to tokenize.
-                # subtle bugs result otherwise. see GH3350
-                token2 = tokenize(name, token if token else var._data, chunks)
-                name2 = f"{name_prefix}{name}-{token2}"
-                return var.chunk(chunks, name=name2, lock=lock)
-            else:
-                return var
-
-        variables = {k: maybe_chunk(k, v, chunks) for k, v in self.variables.items()}
+        variables = {
+            k: _maybe_chunk(k, v, chunks, token, lock, name_prefix)
+            for k, v in self.variables.items()
+        }
         return self._replace(variables)
 
     def _validate_indexers(
@@ -4260,7 +4271,6 @@ class Dataset(Mapping, ImplementsDatasetReduce, DataWithCoords):
         keep_attrs: bool = None,
         keepdims: bool = False,
         numeric_only: bool = False,
-        allow_lazy: bool = None,
         **kwargs: Any,
     ) -> "Dataset":
         """Reduce this dataset by applying `func` along some dimension(s).
@@ -4335,7 +4345,6 @@ class Dataset(Mapping, ImplementsDatasetReduce, DataWithCoords):
                         dim=reduce_dims,
                         keep_attrs=keep_attrs,
                         keepdims=keepdims,
-                        allow_lazy=allow_lazy,
                         **kwargs,
                     )
 
diff --git a/xarray/core/variable.py b/xarray/core/variable.py
index f4ced459..57dfe3d9 100644
--- a/xarray/core/variable.py
+++ b/xarray/core/variable.py
@@ -1582,7 +1582,6 @@ class Variable(
         axis=None,
         keep_attrs=None,
         keepdims=False,
-        allow_lazy=None,
         **kwargs,
     ):
         """Reduce this array by applying `func` along some dimension(s).
@@ -1624,24 +1623,14 @@ class Variable(
         if dim is not None:
             axis = self.get_axis_num(dim)
 
-        if allow_lazy is not None:
-            warnings.warn(
-                "allow_lazy is deprecated and will be removed in version 0.16.0. It is now True by default.",
-                DeprecationWarning,
-            )
-        else:
-            allow_lazy = True
-
-        input_data = self.data if allow_lazy else self.values
-
         with warnings.catch_warnings():
             warnings.filterwarnings(
                 "ignore", r"Mean of empty slice", category=RuntimeWarning
             )
             if axis is not None:
-                data = func(input_data, axis=axis, **kwargs)
+                data = func(self.data, axis=axis, **kwargs)
             else:
-                data = func(input_data, **kwargs)
+                data = func(self.data, **kwargs)
 
         if getattr(data, "shape", ()) == self.shape:
             dims = self.dims


Test patch:
diff --git a/xarray/tests/test_backends.py b/xarray/tests/test_backends.py
index 81664737..b2e2ad3e 100644
--- a/xarray/tests/test_backends.py
+++ b/xarray/tests/test_backends.py
@@ -907,11 +907,13 @@ class CFEncodedBase(DatasetIOBase):
         ve = (ValueError, "string must be length 1 or")
         data = np.random.random((2, 2))
         da = xr.DataArray(data)
-        for name, e in zip([0, (4, 5), True, ""], [te, te, te, ve]):
+        for name, (error, msg) in zip([0, (4, 5), True, ""], [te, te, te, ve]):
             ds = Dataset({name: da})
-            with raises_regex(*e):
+            with pytest.raises(error) as excinfo:
                 with self.roundtrip(ds):
                     pass
+            excinfo.match(msg)
+            excinfo.match(repr(name))
 
     def test_encoding_kwarg(self):
         ds = Dataset({"x": ("y", np.arange(10.0))})
@@ -4307,17 +4309,17 @@ class TestValidateAttrs:
             ds, attrs = new_dataset_and_attrs()
 
             attrs[123] = "test"
-            with raises_regex(TypeError, "Invalid name for attr"):
+            with raises_regex(TypeError, "Invalid name for attr: 123"):
                 ds.to_netcdf("test.nc")
 
             ds, attrs = new_dataset_and_attrs()
             attrs[MiscObject()] = "test"
-            with raises_regex(TypeError, "Invalid name for attr"):
+            with raises_regex(TypeError, "Invalid name for attr: "):
                 ds.to_netcdf("test.nc")
 
             ds, attrs = new_dataset_and_attrs()
             attrs[""] = "test"
-            with raises_regex(ValueError, "Invalid name for attr"):
+            with raises_regex(ValueError, "Invalid name for attr '':"):
                 ds.to_netcdf("test.nc")
 
             # This one should work
@@ -4328,12 +4330,12 @@ class TestValidateAttrs:
 
             ds, attrs = new_dataset_and_attrs()
             attrs["test"] = {"a": 5}
-            with raises_regex(TypeError, "Invalid value for attr"):
+            with raises_regex(TypeError, "Invalid value for attr 'test'"):
                 ds.to_netcdf("test.nc")
 
             ds, attrs = new_dataset_and_attrs()
             attrs["test"] = MiscObject()
-            with raises_regex(TypeError, "Invalid value for attr"):
+            with raises_regex(TypeError, "Invalid value for attr 'test'"):
                 ds.to_netcdf("test.nc")
 
             ds, attrs = new_dataset_and_attrs()
@@ -4668,3 +4670,24 @@ def test_extract_zarr_variable_encoding():
         actual = backends.zarr.extract_zarr_variable_encoding(
             var, raise_on_invalid=True
         )
+
+
+@requires_h5netcdf
+def test_load_single_value_h5netcdf(tmp_path):
+    """Test that numeric single-element vector attributes are handled fine.
+
+    At present (h5netcdf v0.8.1), the h5netcdf exposes single-valued numeric variable
+    attributes as arrays of length 1, as oppesed to scalars for the NetCDF4
+    backend.  This was leading to a ValueError upon loading a single value from
+    a file, see #4471.  Test that loading causes no failure.
+    """
+    ds = xr.Dataset(
+        {
+            "test": xr.DataArray(
+                np.array([0]), dims=("x",), attrs={"scale_factor": 1, "add_offset": 0}
+            )
+        }
+    )
+    ds.to_netcdf(tmp_path / "test.nc")
+    with xr.open_dataset(tmp_path / "test.nc", engine="h5netcdf") as ds2:
+        ds2["test"][0].load()
diff --git a/xarray/tests/test_variable.py b/xarray/tests/test_variable.py
index 08fe0739..cad565f5 100644
--- a/xarray/tests/test_variable.py
+++ b/xarray/tests/test_variable.py
@@ -1570,10 +1570,6 @@ class TestVariable(VariableSubclassobjects):
 
         with raises_regex(ValueError, "cannot supply both"):
             v.mean(dim="x", axis=0)
-        with pytest.warns(DeprecationWarning, match="allow_lazy is deprecated"):
-            v.mean(dim="x", allow_lazy=True)
-        with pytest.warns(DeprecationWarning, match="allow_lazy is deprecated"):
-            v.mean(dim="x", allow_lazy=False)
 
     @pytest.mark.parametrize("skipna", [True, False])
     @pytest.mark.parametrize("q", [0.25, [0.50], [0.25, 0.75]])