From #983 (comment)

With current Shapely, you get an empty tuple for the bounds of an empty geometry:

>>> from shapely.geometry import Polygon
>>> p = Polygon()
>>> p.bounds
()

In PyGEOS, we return 4 NaNs, instead:

In [1]: import pygeos

In [2]: pygeos.bounds(pygeos.Geometry("POLYGON EMPTY"))
Out[2]: array([nan, nan, nan, nan])

For the vectorized ufuncs, I think this makes sense. And actually also the only possible solution, because when applying this operation on multiple points, we need to be able to put the results in a single array with consistent dimension:

In [3]: pygeos.bounds([pygeos.box(0, 0, 1, 1), pygeos.Geometry("POLYGON EMPTY")])
Out[3]: 
array([[ 0.,  0.,  1.,  1.],
       [nan, nan, nan, nan]])

So the question is then: for the scalar .bounds attribute on the geometry objects, do we keep the current empty tuple () result, or do we chagne to 4-length tuple (nan, nan, nan, nan) as well?

I would personally propose to change to 4 NaNs for the scalar attribute: 1) consistency with the vectorized operation, 2) consistency in the length of the returned value (eg p.bounds[0] will always work)