[RFC] Make .X a layer

[gtca]

This is an implementation proposal draft for unifying .X and layers as mentioned e.g. in #244. Expanding on previous issues, this one tracks concrete implementation details for achieving the end goal.

At this stage, feedback on the topic is appreciated — either here or in #244.

.X attribute as a reference

The proposal is to use adata.X as a reference to one of the layers. In the text of the proposal below, *.X or *adata.X will denote the data in the layer that adata.X points to (akin to pointer dereferencing).

Benefits

Some benefits include:

• Reduce the implementation complexity as .X and .layers items store conceptually identical entities.
• Reduce memory footprint in common workflows. E.g. when having layers with spliced and unspliced counts, adata.X will just point to one of them instead of containing its copy.
• Make it easier to implement idempotent APIs (to be expanded in other issues).

Options

1. Explicit layer name

.X links to an explicit layer. This reference can be changed by a function that is run on an AnnData object. E.g. a normalisation function can create a new layer and update the reference instead of modifying the values in-place. In-place APIs will continue to modify the values in the *.X.

When *.X is deleted, .X can be either None or, less favourably, point to the layer that was defined previously in case there's a memory of which layers were linked to in .X before.

2. Last layer

.X can link to adata.layers[list(adata.layers.keys())[-1]]. This is easier to reason about but might be an issue when there's no explicit order in layers (e.g. raw counts -> normalised counts -> scaled counts) and new layers are added (e.g. unspliced counts).

Additional methods

There should be a way to explicitly change which layer .X points to.

Alternatives

An alternative, as mentioned in #244, would be to have .X as a specific layer. While this is simpler in the short-term, that might not provide some of the benefits listed above.

Effect for the end user

Taking the current scanpy API as an example where the matrix is modified in-place, nothing changes for the end user:

sc.pp.scale(adata)
# modifies values stored in the layer
# that .X links to

On disk

For HDF5, a non-breaking change should be achievable with the HDF5 references functionality. It is to be checked that R (rhdf5, hdf5r) and Julia HDF5 libraries are able to work with references.

For Zarr, feedback is needed if this is possible. It is probably not an issue if it's not as Zarr storage for AnnData is less mature/common.

Potential

This is to be expanded in other issues to come.

Idempotent analysis API

In practice, the scaling example above typically looks like this together with its context:

adata.layers["lognorm"] = adata.X.copy()
sc.pp.scale(adata)

This feature would allow to implement idempotent APIs, i.e.

sc.pp.scale(adata, in_layer="lognorm", out_layer="scaled")
# .X now links to .layers["scaled"]

that would not modify counts, and, importantly, will not affect the end user API provided that the appropriate default arguments for in_layer and out_layer are chosen.

[grst]

Probably makes sense to implement this before we move to supporting Awkward Arrays (or fwiw other data types) in adata.X (see also #647).

I was wondering if it would make sense that adata.var and adata.obs become pointers to a "layer" in adata.varm and adata.obsm, too?

[flying-sheep]

I think this is a nicer idea than the one in #244, so we can just follow this proposal instead of mine!

Duplicate of #244