Revise how the domain of an array is specified

[jbms]

This adds support for dimension names (#73) and non-zero origins (#122).

[alimanfoo]

Just to say thanks @jbms.

Reading through #122 is very enlightening. Do we think there are any outstanding objections that would need further discussion before merging here?

[jbms]

I think the main objection was from zarr.jl, in particular I think there is some disagreement as to how interoperability between languages should work as far as index vectors. The zarr.jl model is that (i, j, k) in zarr-python should correspond to (k+1, j+1, I+1) in zarr.jl., while I would prefer that zarr.jl users just define an array with a Fortran order layout and inclusive_min of (1, 1, 1) if that is what is desired, so that there is no confusion or inconsistency with indexing. Still I think there are some options for how this feature could be supported in zarr.jl while retaining the option to do that coordinate transform.

[jbms]

Perhaps since no objections have been raised in the last few days it makes sense to merge this and address any objections as part of the general ZEP 1 discussion?

I have some additional changes I'd like to propose that build on this.

[meggart]

I have just read the proposed changes and as you already mentioned I would be quite unhappy when we force offset indexing into the core data model. I tried to read #122 carefully and still do not understand why these types of transformations have to be part of the core data model and can not be defined in an extension.

If I implement the deafults you suggest in Julia, imagine I save and publish a dataset with Fortran order and define a domain using 1-based indexing. Now everyone loading my dataset from python/Javascript/C++ etc will have to access my data using 1-based indexing, just because this is what I used in my Julia code? And in addition they have to remember to revert their loop-dimensions for cache-friendly iteration through the dataset.

And the same applies for the other way: All Julia users will have to deal with 0-based OffsetArrays and PermutedDimsArrays just to make sure to have a 1-1 mapping between array domains between different programming languages.

I am in full support of the idea to allow negative chunk indices and and for allowing custom offsets, I even use an _ARRAY_OFFSET attribute for the Zarr arrays I am writing in YAXArrays.jl, but I would be worried if we enforce a particular indexing semantics onto all existing implementations.

[rabernat]

I am also opposed to adding offset indexing into the core Zarr data model, for the reasons I outlined in #122 (comment).

I recognize that this feature is very important to some user communities. We should definitely find a way to to address it. But I don't think core specification is the right place.

The core specification describes the minimum feature set that all implementations must support. Under the core specification, I believe that each implementation should expose indexing on Zarr arrays using its own natural language-specific indexing convention, i.e. 0-based C-order for python / C, 1-based F-order for Julia / Fortran, etc.

Instead, an extension could be defined which defines a transformation between the raw index space and the offset index space. Implementations could then choose to expose the offset indexes or the raw indexes.

@jbms - can you explain why you think this feature needs to be in the core spec as opposed to an extension?

[jbms]

Thanks for your feedback @rabernat and @meggart.

To be clear, this PR only adds support for non-zero origins, but does not add support for arbitrary dimension orders. I was planning to create a PR for arbitrary dimension order support, but I think it would be simpler if for now we just discuss non-zero origin support.

As I see it, being able to consistently and unambiguously refer to positions or regions of an array is a core component of an array storage standard. While I can understand that this is not important for some applications, I think that this tends to be more important as arrays become too large to fit in memory, which is precisely the class of applications for which a chunked array format like zarr is most useful.

@rabernat asks why non-zero origins support isn't better supported as an extension. I see several possibilities:

• If it is added as a must_understand = false extension, and some implementations don't support it, then rather than it helping to clarify the bounds of the data, it may likely lead to more user confusion, as the indexing behavior may vary across implementations depending on whether the extension is supported. A must_understand = false extension also could not support dimensions where the grid_origin is greater than the inclusive_min, which is useful to allow an array to expand "downwards".

• If it is added as a must_understand = true extension, and some/most implementations don't support it, then we lose interoperability, and therefore lose much of the benefit of using zarr at all.

• If it is added as a must_understand = true extension and all/most implementations support it, then it is effectively no different from having it as part of the core specification.

While this feature certainly has implications on the zarr data model, I don't think it imposes a significant implementation burden --- it is just two additional optional parameters that must be managed, and taken into account in indexing operations.

@meggart raises the concern that this feature might require support for in-memory arrays with non-zero origins in every zarr implementation, e.g. OffsetArray in Julia. But I think that is not necessarily the case. For example, TensorStore supports non-zero origins but in its Python API, reading from a TensorStore simply returns a normal NumPy array (corresponding to whichever region was requested). Similarly, I would imagine that reading a portion of a Zarr v3 array in Zarr.jl could just return a regular Julia array, which could then be processed in the normal way. In general while in some ways a zarr array and an in-memory array are similar and could be used interchangeably, they differ in many important ways (latency, efficient access patterns, fallibility of operations). Therefore it seems reasonable that they may not have identical APIs.

[rabernat]

I'm curious to hear the Unidata perspective. @WardF & @DennisHeimbigner: would would be the implication of putting this feature in the core specification for NetCDF-C? How would you deal with the fact that Zarr arrays could potentially have non-zero origin as in this proposal?

[DennisHeimbigner]

We might be able to at least pass info to the user by exposing an attribute that
gives the offsets from the zarr file, even if we do not use them.

[DennisHeimbigner]

I suppose that we can handle this similar to the existing fortran solution.
I presume that the transform is only applied to the array coordinates
passed in thru the API. So, the underlying storage layout is unaffected.
So given an array(d1,d2,d3) and offset of (x,y,z), and request for data array(p1,p2,p3)
we transform to a request for array(p1-x,p2-y,p3-z) and all else proceeds internally
as normal. Correct?

[jbms]

I suppose that we can handle this similar to the existing fortran solution. I presume that the transform is only applied to the array coordinates passed in thru the API. So, the underlying storage layout is unaffected. So given an array(d1,d2,d3) and offset of (x,y,z), and request for data array(p1,p2,p3) we transform to a request for array(p1-x,p2-y,p3-z) and all else proceeds internally as normal. Correct?

That is correct, but note that the relevant offset (x, y, z) is the grid_origin, not the inclusive_min value for each dimension. The inclusive_min is only used for bounds checking, and does not affect how the data is retrieved.

[rabernat]

Summarizing some of the discussion that occurred about this on yesterday's call:

I argued that the Zarr specification essentially remain agnostic about indexing. Different languages have different index idioms for addressing the elements in ND-arrays. The Zarr implementation in a specific language should use that language's conventions to access these elements. So the element I address as foo[0, 1] in Python would be foo[2, 1] in Julia, reflecting both the different origin and the different axis order (see #126). It is important to recognize that these are the same data structures as described by different languages, just like an apple and a manzana are the same fruit in English and Spanish, despite different spelling. Doing otherwise would hurt adoption and interoperability of Zarr within that language community by exposing arrays with unconventional / unexpected behavior. A Zarr array should look and feel as much as possible like a normal array (whatever that means) from that language.

If we accept that the spec should remain agnostic about indexing conventions, it makes it very hard to write the spec, because at the end of the day, the spec must refer to specific elements in its examples. To mitigate this, we simply state that the spec adopts zero-based C-style indexing when describing specific array elements, and that other implementations may use other indexing conventions.

I argued that anything other than this raw, idiomatic indexing must be considered a coordinate transformation from some coordinate space to raw index space. We can denote this N-dimensional raw index space as

$$ \mathbb{I}^N := \{ (i_1, i_2, ..., i_N) | i_n \in \mathbb{Z}^{0+} \} $$

(using of course the 0-origin convention.)

All other coordinate system are ultimately mappings from some other space $X$ to index space: $f : X \to \mathbb{I}^N$. In this notation, the current PR proposes a mapping

$$ f : \mathbb{Z}^N \to \mathbb{I}^N $$

where the function $f$ is just the addition of an offset, i.e. $f(i) = i + C$.

The affine transforms proposed above could obviously be framed in the same way. CF Conventions on coordinates describe a whole system for mapping discrete coordinate values (stored in other arrays within a group) to raw index space. Xarray (and many other similar software applications) leverage this convention to support label based indexing.

It is obviously my opinion that we should not bring any coordinate transformations into the core Zarr spec. In fact, I don't even think we need extensions for them. As the CF conventions show, such coordinate transforms can be described completely by appropriate attributes, without any spec changes needed. For non-zero origin, it could be as simple as an array attribute like:

index_origin: (-10000, 3000, 0)

If we accept this, the challenges (for Jeremy and others who want this feature) are:

• Where should the convention for non-zero origin be written down? (I would be happy to see Zarr host such attribute conventions, separate from the spec itself.)
• What layer of the stack is responsible for encoding / decoding these conventions, if not Zarr itself? I recognize this is a serious problem for communities who work with "raw" zarr directly, without an Xarray-like intermediary. Julia, for example, has OffsetArrays.jl/. We would need a package like this for Python, and a dispatch mechanism from Zarr.

One place where we have a vaguely similar situation is with fill_value. Zarr does specify a fill_value in the core spec. NetCDF / CF conventions have a separate fill_value attribute, which we can put in Zarr attributes. Just today we were working on a bug that arose from the fact that there are two places in the stack that address fill values, and these can become inconsistent. A very clean separation between the requirements of the storage layer (Zarr / HDF / etc) and the higher-level data model (NetCDF / Xarray / etc) helps avoid these problems. That's the basic reason why I am opposed to adding more higher-level features related to coordinates into the core Zarr spec.

[jbms]

I can see that there is a strong opposition to non-zero origins in zarr, so I won't continue pushing it.

As suggested, I can indeed implement support outside of the spec in Neuroglancer and TensorStore using user-defined attributes, and will try to coordinate with other implementations that are interested, such as WebKnossos.

The one aspect of this proposal that is not as conveniently implemented outside of the core spec is support for growing a dimension in the negative direction, or similarly, allowing dimensions to be unbounded in the negative direction. (To use the mathematical notation, the most significant change in this PR is that Zarr arrays would be indexed by $Z^N$ rather than $I^N$. Since negative chunks are allowed, the grid_origin parameter simply maps from $Z^N$ to $Z^N$.)

However, that can still be accomplished by choosing an initial origin that is some large number, e.g. $2^{62}$. The downside is that it would be even more awkward to access the array with zarr implementations that are unaware of the origin metadata, but that may be acceptable when an unbounded domain/resizing in the negative direction is required.

There is one remaining portion of this PR that may be less controversial, and which hasn't been discussed yet: dimension names. This is a feature that is very easy to implement outside of the core spec, e.g. a simple "dimension_names" user-defined attribute, but there are some advantages to making it part of the core spec:

• Dimension labels are widely used by various tools and specifications that build on top of zarr, e.g. xarray, OME, nczarr, tensorstore, neuroglancer.
• If some implementations such as Zarr.jl choose to reverse the dimension order, dimension labels are particularly important to avoid ambiguity.
• Implementation effort is extremely low, in fact implementations that don't wish to support it could ignore the labels altogether without causing any problems.
• If not specified in the core spec, there is a risk that multiple incompatible conventions may be developed. In particular, there is the risk that users of implementations that reverse the dimensions may end up storing their dimension labels in reverse order.

[rabernat]

The one aspect of this proposal that is not as conveniently implemented outside of the core spec is support for growing a dimension in the negative direction, or similarly, allowing dimensions to be unbounded in the negative direction.

I agree--this is a tough one. My main objection was that we don't want to mix the notion of coordinates into the core spec. But I do agree that it would be great to allow extending arrays in both directions. I think that can be done without introducing coordinates. Maybe we could support this via an optional extension that would operate at the chunk level, i.e. creating chunks like -1.0.0. Have you though about how negative indexes would work from within python at an API level?

---

I am 100% in favor of dimension labels.

• Implementation effort is extremely low, in fact implementations that don't wish to support it could ignore the labels altogether without causing any problems.

Doesn't this kind of imply it should be an extension?

If not specified in the core spec, there is a risk that multiple incompatible conventions may be developed.

In fact, this has already happened with xarray and nczarr. But I'm convinced that could have been avoided also with a clearly documented extension.

I think the only question for me is--core spec vs. extension? You seem to be suggesting that things that are not in the core spec will just be ignored. I was hoping that extensions would be a little more binding and widely adopted. Named dimensions has always been held up as an example of what an extension would do. Very curious to get thoughts from others here.

[d-v-b]

I agree that these augmentations to array semantics are both extremely useful and should be part of an extension. If everyone deems these augmentations essential and we observe that essentially every zarr array in the wild uses the extension, then it can be later added to the core spec.

[jbms]

The one aspect of this proposal that is not as conveniently implemented outside of the core spec is support for growing a dimension in the negative direction, or similarly, allowing dimensions to be unbounded in the negative direction.

I agree--this is a tough one. My main objection was that we don't want to mix the notion of coordinates into the core spec. But I do agree that it would be great to allow extending arrays in both directions. I think that can be done without introducing coordinates. Maybe we could support this via an optional extension that would operate at the chunk level, i.e. creating chunks like -1.0.0. Have you though about how negative indexes would work from within python at an API level?

I think there are a few different options:

a. Don't support the "negative index means counting from the end" convention at all.
b. Support the "negative index means counting from the end" convention by default, but don't support it on arrays for which negative chunks/indices are supported.
c. Use a separate accessor, e.g. arr.nindex[-1, -5:3] to disable the "negative index means counting from the end" convention.
d. Require a special option when opening the array, e.g. zarr.open(..., allow_negative_indices=True).

The TensorStore Python API uses option (a). For zarr-python I think option (b) might be the most practical.

Presumably with such an extension, the shape attribute in the metadata would be replaced by inclusive_min and exclusive_max. (Alternatively shape could still be used along with inclusive_min, except that it is then problematic to have an unbounded dimension.) Then the only functionality difference between that extension and this PR would be that this PR also includes a grid_origin option.

Do I understand correctly that you are supportive of adding support for such an extension to zarr-python?

I am 100% in favor of dimension labels.

• Implementation effort is extremely low, in fact implementations that don't wish to support it could ignore the labels altogether without causing any problems.

Doesn't this kind of imply it should be an extension?

I think it would be fine as an extension as well --- I don't think there would be much practical difference.

I can can create a PR that adds it as an extension if that is preferable. The main thing would be to standardize it early so that multiple incompatible ways to specify dimension labels aren't developed.

I imagined that the core spec could also list some features as optional, but if the idea is that nothing in the core spec is optional then I can see an extension might be a better fit. Though it also seems plausible that every implementation would support at least the ability to query the list of dimension labels of an open array, and set the dimension labels when creating an array.

If not specified in the core spec, there is a risk that multiple incompatible conventions may be developed.

In fact, this has already happened with xarray and nczarr. But I'm convinced that could have been avoided also with a clearly documented extension.

I think the only question for me is--core spec vs. extension? You seem to be suggesting that things that are not in the core spec will just be ignored. I was hoping that extensions would be a little more binding and widely adopted. Named dimensions has always been held up as an example of what an extension would do. Very curious to get thoughts from others here.

I think that after the initial version of the spec is released, extensions become an important technical measure for adding functionality in a way that explicitly indicates whether old implementations unaware of the extension should either ignore it or fail immediately.

Prior to the release of the initial version of the spec, for optional features that don't otherwise require e.g. a data type/codec/chunk layout/storage transformer identifier, it is somewhat arbitrary whether it is an extension.

I agree that these augmentations to array semantics are both extremely useful and should be part of an extension. If everyone deems these augmentations essential and we observe that essentially every zarr array in the wild uses the extension, then it can be later added to the core spec.

I am perfectly happy for this (negative indices) functionality to be an extension rather than in the core spec; but I think it would be significantly more useful if support for it would also be accepted for inclusion in zarr-python.

[d-v-b]

I don't understand what the trouble is for growing the array in the "negative" direction. Why can't zarr just copy the numpy pad API? I see no need for breaking the conventions for negative indexing here (just as np.pad doesn't break any array indexing conventions).

[joshmoore]

Why can't zarr just copy the numpy pad API?

This would require moving every chunk file, no? And on s3 that would be a copy+delete.

[jbms]

@d-v-b raises an interesting point --- if grid_origin is supported (but the lower bound is still always 0) then it would indeed be possible to support an API like pad. Then I could still support negative indices through some additional extension for which zarr-python support might not be acceptable. This approach may raise some issues with concurrent access, though, depending on the use case.

[d-v-b]

This would require moving every chunk file, no? And on s3 that would be a copy+delete.

Is it baked into the spec that the names of chunks in storage must be positive integers? The name of a chunk in storage is different from how it is addressed from zarr-python, so you could allow chunks with negative integers in their names (provided the storage backend allows it), but the first chunk is always addressed by array[0, 0, ...].

However, as noted earlier, this would require explicitly tracking the name of the first chunk in some metadata.

[rabernat]

I imagined that the core spec could also list some features as optional, but if the idea is that nothing in the core spec is optional then I can see an extension might be a better fit.

☝️ this strikes me as a really core issue. Curious to get @alimanfoo's take on it.

[jstriebel]

Since this discussion stalled, I've tried to summarize the different points here, added my and @normanrz's comments and added some proposals how to proceed, to keep progressing with this PR. I hope I didn't oversee anything in the discussion above and am sorry for the long text 🙈

Summary

IMO there are six different points adressed in the discussion (and PR):

1. Adding (optional) names to dimensions. There seems to be agreement that his is helpful. More details can also be found in Dimension names as core array metadata #73.
2. Adding min and max per dimension, instead of starting at 0 and having a shape (Support for non-zero origin #122). There seem to be different opinions if this should be part of the core model, or if arrays should always start at 0 and can possibly be transformed additionally. The following data models are thinkable (only considering zarr indexing here, not any higher-level libs):
   1. All data should start at 0, shape specifies the maximum valid index + 1.
   2. All data can start at any positive coordine, shape specifies the maximum valid index + 1.
      In this case the minimum is unknown.
   3. Same as 2, but define the minimum as well.
   4. Same as 3, but also allow negative coordinates.
      In all of those cases there may be additional transformations in a fixed metadata field, that either must be applied for indexingt (see 3) or that can be applied by higher-level libraries and viewers (see 4.)
3. An indexing offset could change how data is stored in comparison how it is indexed (e.g. find all indices X where we'd normally store position X+offset). Comment: For an arbitrary offset this would need to be part of the core spec. When only allowing chunk-wise offsets (so offset must be divisible by the chunk-size), this could be handled by a storage transformer.
4. Domain transformations can be added to fit the data into a domain-space which might differ from the addressing space. zarr would not directly use this for adressing. It might become a metadata convention (see 6.).
5. It would be great if left-padding data is possible efficiently, meaning without renaming all written chunks. This is depending on 2., as some of those don't allow to insert chunks to the "left" (smaller coordinates that the given minimum):
   • 2.i. This does not allow efficient left-padding.
   • 2.ii. and iii. This does allow efficient left-padding as long as there is enough space left between 0 and the actual data.
   • 2.iv. This does allow efficient left-padding in all cases.
     (For this discussion I'll ignore changes to the leftmost chunk, this is assuming that more is added than what fit's into the leftmost chunks)
     One argument against negative indices (in the API) was, that
     Another option is to use the indexing offset from 3. to remap the chunks, this allows to pad further to the left than what previously was 0. also for the cases 2.i.-2.iii. In the storage chunk-indices might be read from negative numbers, but in the API they would appear positive. Inserting data and changing this offset would then effectively move all the previous data to the "right", indexed coordinates would change.
6. It might be useful to standardize where to capture metadata that might not be used by zarr itself for indexing/reading/writing as metadata/attributes conventions. This seems to be different from extensions, as they don't change behavior, but are recommendations to support compatibility between higher-level libs and applications. 1 and 4 seem to be such candidates.

Comments

on min/max (2.)

I'll try to summarize @normanrz's and my thoughts on the proposal to add a minimum:
Currently (as of v2), data starting at 0 (2.i.) seems to be expected, but actually starting later (2.ii. is possible), and also needed to efficiently allow left-padding (5.). Another use case for not starting the data at 0 is when having multiple arrays in a shared space, but with different bounding boxes, and possibly offsets between them. It would be favourable, that they have shared indexing coordinates a) for convenience, and b) to ensure the same chunking granularity. An example: 2D arrays X and Y are both in the same domain (e.g. microscopy data for X and a segmentation of a subset for Y), both with chunk-size 32, 32. Let's assume X starts at 0, 0, but Y at 42, 42. It would be great if a) we could address the same array areas with the same zarr indices, b) the same chunk-borders are used and c) we could know the valid range of Y. For such a use-case, saving the minimum index as well would be great. Defining an indexing offset for Y of 42,42 would only allow to store the chunks in a coherent space, but the indexing via zarr would still differ between X and Y.

In general, it seems that saving the minimum as well is similar to saving the shape. shape is mostly used to know what part of an array is expected to be defined, and to allow to read "the whole array" in the API. How to use the minimum or maximum/shape is up the implementation. E.g. in the case of zarr-python the API is modelled to be close to numpy-array atm, so indexing beyond the shape implicitly returns the data until the shape (zarr.ones((10, 10))[:20, :20].shape and np.ones((10, 10))[:20, :20].shape both are (10, 10)). Another implementation might work differently, e.g. returning the fallback value beyond the array. I think this behavior does not need to be specified in the spec, rather explicitly stated that everything beyond the array border is undefined. Even when specifying the minimum in the metadata, zarr-pyhon might still decide to return the array starting from 0 in the numpy-compatible api, but might also expose the minimum. To rephrase it: the index-minimum would just be another metadata attribute the implementation can decide to use to present an API as expected.

Regarding negative indices: Since negative indices would create problems in APIs, I'd rather avoid them. Left-padding beyond 0 can still be achieved by an indexing offset (4.).

In short: I'd just change the allowed ranges for min/max to be positive, otherwise the changes LGTM.

on domain names (1.), domain transformations (3.) and metadata/attribute conventions (6.)

The only question left for this PR IMO is if this should be part of the attributes or the general metadata. I'd tend to move everything that is not used for indexing/reading/writing into attributes, but domain names might be important enough to be an exception.

Domain transformations don't seem to be well-defined atm and don't change anything in the core zarr mechanics, so I'd rather defer this until multiple domains converge on similar terms (e.g. compare ome/ngff#94). I think such conventions can easily be added later.

Next steps

I think we should separate the different points a bit and concentrate on the actual PR proposal at hand. Going through the different points:

1. dimension names Let's discuss this further in this PR.
2. min and max As well. It might make sense to separate this from 1. if there is no agreement.
3. indexing offset This would be possible via the storage transformer extension point, when only allowing chunk-divisible offsets. I'd suggest to use a separate issue if someone wants to propose this, or propose arbitrary offsets in the v3 core. We should clarify if the arbitrary offset would be an alternative to the index-minimum or if this might be an addition.
4. domain transformations This seems to be a separate topic which is probably better off in a separate issue.
5. left-padding Might be worth a separate issue, as this might not be fully resolved by this PR. It might even get it's own section in the spec (e.g. under implementation), as it's quite a complicated detail.
6. metadata/attributes conventions: I'll open a separate issue for those.

[jbms]

I've pulled out just the addition of dimension names into a separate PR:

#162

[jstriebel]

@jbms and I decided to close this since #162 is in place for dimension names and non-zero origins (#122) is proposed as an extension, not being part of v3 core