Repository: single-cell-data/TileDB-SOMA (language: python)

Issue:
[c++/python/r] Revert C++20 to C++17 on the `release-1.15` branch
Revert C++20 changes to unblock release


Discussion / hints:
## [Codecov](https://app.codecov.io/gh/single-cell-data/TileDB-SOMA/pull/3389?dropdown=coverage&src=pr&el=h1&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=single-cell-data) Report
All modified and coverable lines are covered by tests :white_check_mark:
> Project coverage is 85.96%. Comparing base [(`684cd69`)](https://app.codecov.io/gh/single-cell-data/TileDB-SOMA/commit/684cd69d83547a78db8d44b4f51a912c2bc42506?dropdown=coverage&el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=single-cell-data) to head [(`793faea`)](https://app.codecov.io/gh/single-cell-data/TileDB-SOMA/commit/793faeaa8a51326c55b693f4adf5238b5415d76c?dropdown=coverage&el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=single-cell-data).
> Report is 1 commits behind head on release-1.15.

<details><summary>Additional details and impacted files</summary>


```diff
@@               Coverage Diff                @@
##           release-1.15    #3389      +/-   ##
================================================
+ Coverage         85.87%   85.96%   +0.09%     
================================================
  Files                55       55              
  Lines              6143     6143              
================================================
+ Hits               5275     5281       +6     
+ Misses              868      862       -6     
```

| [Flag](https://app.codecov.io/gh/single-cell-data/TileDB-SOMA/pull/3389/flags?src=pr&el=flags&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=single-cell-data) | Coverage Δ | |
|---|---|---|
| [python](https://app.codecov.io/gh/single-cell-data/TileDB-SOMA/pull/3389/flags?src=pr&el=flag&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=single-cell-data) | `85.96% <ø> (+0.09%)` | :arrow_up: |

Flags with carried forward coverage won't be shown. [Click here](https://docs.codecov.io/docs/carryforward-flags?utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=single-cell-data#carryforward-flags-in-the-pull-request-comment) to find out more.

| [Components](https://app.codecov.io/gh/single-cell-data/TileDB-SOMA/pull/3389/components?src=pr&el=components&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=single-cell-data) | Coverage Δ | |
|---|---|---|
| [python_api](https://app.codecov.io/gh/single-cell-data/TileDB-SOMA/pull/3389/components?src=pr&el=component&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=single-cell-data) | `85.96% <ø> (+0.09%)` | :arrow_up: |
| [libtiledbsoma](https://app.codecov.io/gh/single-cell-data/TileDB-SOMA/pull/3389/components?src=pr&el=component&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=single-cell-data) | `∅ <ø> (∅)` | |

</details>

Fix patch:
diff --git a/.github/workflows/python-ci-packaging.yml b/.github/workflows/python-ci-packaging.yml
index 05eda86167..64cff69e38 100644
--- a/.github/workflows/python-ci-packaging.yml
+++ b/.github/workflows/python-ci-packaging.yml
@@ -76,7 +76,7 @@ jobs:
         run: |
           mkdir -p external
           # Please do not edit manually -- let scripts/update-tiledb-version.py update this
-          wget --quiet https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0-rc3/tiledb-linux-x86_64-2.27.0-rc3-8d581f2.tar.gz
+          wget --quiet https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0.rc5/tiledb-linux-x86_64-2.27.0.rc5-2862c30.tar.gz
           tar -C external -xzf tiledb-linux-x86_64-*.tar.gz
           ls external/lib/
           echo "LD_LIBRARY_PATH=$(pwd)/external/lib" >> $GITHUB_ENV
@@ -178,10 +178,10 @@ jobs:
           mkdir -p external
           # Please do not edit manually -- let scripts/update-tiledb-version.py update this
           if [ `uname -m` == "arm64" ]; then
-            wget --quiet https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0-rc3/tiledb-macos-arm64-2.27.0-rc3-8d581f2.tar.gz
+            wget --quiet https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0.rc5/tiledb-macos-arm64-2.27.0.rc5-2862c30.tar.gz
             tar -C external -xzf tiledb-macos-arm64-*.tar.gz
           else
-            wget --quiet https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0-rc3/tiledb-macos-x86_64-2.27.0-rc3-8d581f2.tar.gz
+            wget --quiet https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0.rc5/tiledb-macos-x86_64-2.27.0.rc5-2862c30.tar.gz
             tar -C external -xzf tiledb-macos-x86_64-*.tar.gz
           fi
           ls external/lib/
@@ -274,14 +274,14 @@ jobs:
           if [ `uname -s` == "Darwin" ]; then
             if [ `uname -m` == "arm64" ]; then
                 # Please do not edit manually -- let scripts/update-tiledb-version.py update this
-                wget --quiet https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0-rc3/tiledb-macos-arm64-2.27.0-rc3-8d581f2.tar.gz
+                wget --quiet https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0.rc5/tiledb-macos-arm64-2.27.0.rc5-2862c30.tar.gz
             else
                 # Please do not edit manually -- let scripts/update-tiledb-version.py update this
-                wget --quiet https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0-rc3/tiledb-macos-x86_64-2.27.0-rc3-8d581f2.tar.gz
+                wget --quiet https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0.rc5/tiledb-macos-x86_64-2.27.0.rc5-2862c30.tar.gz
             fi
           else
             # Please do not edit manually -- let scripts/update-tiledb-version.py update this
-            wget --quiet https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0-rc3/tiledb-linux-x86_64-2.27.0-rc3-8d581f2.tar.gz
+            wget --quiet https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0.rc5/tiledb-linux-x86_64-2.27.0.rc5-2862c30.tar.gz
           fi
           tar -C external -xzf tiledb-*.tar.gz
           ls external/lib/
diff --git a/.github/workflows/r-ci.yml b/.github/workflows/r-ci.yml
index 20f98e4aa0..5122f02de0 100644
--- a/.github/workflows/r-ci.yml
+++ b/.github/workflows/r-ci.yml
@@ -99,7 +99,7 @@ jobs:
       #   uses: lukka/get-cmake@latest
 
       #- name: MkVars
-      #  run: mkdir ~/.R && echo "CXX17FLAGS=-Wno-deprecated-declarations -Wno-deprecated" > ~/.R/Makevars
+      #  run: mkdir ~/.R && echo "CXX20FLAGS=-Wno-deprecated-declarations -Wno-deprecated" > ~/.R/Makevars
 
       #- name: Build and install libtiledbsoma
       #  run: sudo scripts/bld --prefix=/usr/local
diff --git a/apis/python/notebooks/README.md b/apis/python/notebooks/README.md
new file mode 100644
index 0000000000..fc16300e21
--- /dev/null
+++ b/apis/python/notebooks/README.md
@@ -0,0 +1,1 @@
+If you add notebooks here, please also update the symlinks in `/doc/source/notebooks`.
diff --git a/apis/python/notebooks/tutorial_soma_shape.ipynb b/apis/python/notebooks/tutorial_soma_shape.ipynb
index 80c57c3e5d..20c0e94d95 100644
--- a/apis/python/notebooks/tutorial_soma_shape.ipynb
+++ b/apis/python/notebooks/tutorial_soma_shape.ipynb
@@ -11,7 +11,9 @@
     "\n",
     "In this notebook, we'll go through how you use shapes for the dataframes and arrays within your SOMA experiments, when and how you can resize, and options for experiments created before TileDB-SOMA 1.15.\n",
     "\n",
-    "The dataset used is from Peripheral Blood Mononuclear Cells (PBMC), which is freely available from 10X Genomics. "
+    "The dataset used is from Peripheral Blood Mononuclear Cells (PBMC), which is freely available from 10X Genomics.\n",
+    "\n",
+    "(Please also see the [Academy tutorial](https://cloud.tiledb.com/academy/structure/life-sciences/single-cell/tutorials/shapes/).)"
    ]
   },
   {
diff --git a/apis/python/requirements_spatial.txt b/apis/python/requirements_spatial.txt
index d91e347662..398dd9fef8 100644
--- a/apis/python/requirements_spatial.txt
+++ b/apis/python/requirements_spatial.txt
@@ -3,4 +3,4 @@ tifffile
 pillow
 spatialdata>=0.2.5
 xarray
-dask
+dask<=2024.11.2
diff --git a/apis/python/src/tiledbsoma/_common_nd_array.py b/apis/python/src/tiledbsoma/_common_nd_array.py
index 885cc8f136..9b28ce5b3c 100644
--- a/apis/python/src/tiledbsoma/_common_nd_array.py
+++ b/apis/python/src/tiledbsoma/_common_nd_array.py
@@ -13,7 +13,7 @@
 from typing_extensions import Self
 
 from ._soma_array import SOMAArray
-from ._types import OpenTimestamp
+from ._types import OpenTimestamp, StatusAndReason
 from .options._soma_tiledb_context import (
     SOMATileDBContext,
 )
@@ -51,9 +51,9 @@ def create(
                 None)``, as the sequence length determines the number of dimensions
                 N in the N-dimensional array.
 
-                For :class:`SparseNDArray` only, if a slot is None, then the maximum
-                possible int32 will be used.  This makes a :class:`SparseNDArray`
-                growable.
+                For :class:`SparseNDArray` only, if a slot is None, then the minimum
+                possible range will be used.  This makes a :class:`SparseNDArray`
+                growable using ``resize``.
             platform_config:
                 Platform-specific options used to create this array.
                 This may be provided as settings in a dictionary, with options
@@ -84,6 +84,39 @@ def create(
         """
         raise NotImplementedError("must be implemented by child class.")
 
+    def resize(
+        self, newshape: Sequence[Union[int, None]], check_only: bool = False
+    ) -> StatusAndReason:
+        """Increases the shape of the array as specfied. Raises an error if the new
+        shape is less than the current shape in any dimension. Raises an error if
+        the new shape exceeds maxshape in any dimension. Raises an error if the
+        array doesn't already have a shape: in that case please call
+        tiledbsoma_upgrade_shape. If ``check_only`` is ``True``, returns
+        whether the operation would succeed if attempted, and a reason why it
+        would not.
+
+        Lifecycle:
+            Maturing.
+        """
+        if check_only:
+            return self._handle.tiledbsoma_can_resize(newshape)
+        else:
+            self._handle.resize(newshape)
+            return (True, "")
+
+    def tiledbsoma_upgrade_shape(
+        self, newshape: Sequence[Union[int, None]], check_only: bool = False
+    ) -> StatusAndReason:
+        """Allows the array to have a resizeable shape as described in the TileDB-SOMA
+        1.15 release notes.  Raises an error if the new shape exceeds maxshape in
+        any dimension. Raises an error if the array already has a shape.
+        """
+        if check_only:
+            return self._handle.tiledbsoma_can_upgrade_shape(newshape)
+        else:
+            self._handle.tiledbsoma_upgrade_shape(newshape)
+            return (True, "")
+
     @property
     def shape(self) -> Tuple[int, ...]:
         """Returns capacity of each dimension, always a list of length ``ndim``.
diff --git a/apis/python/src/tiledbsoma/_dataframe.py b/apis/python/src/tiledbsoma/_dataframe.py
index ec6b02c0ee..973b3e01d7 100644
--- a/apis/python/src/tiledbsoma/_dataframe.py
+++ b/apis/python/src/tiledbsoma/_dataframe.py
@@ -28,7 +28,6 @@
 from . import pytiledbsoma as clib
 from ._constants import SOMA_JOINID
 from ._exception import SOMAError, map_exception_for_create
-from ._query_condition import QueryCondition
 from ._read_iters import TableReadIter
 from ._soma_array import SOMAArray
 from ._tdb_handles import DataFrameWrapper
@@ -172,16 +171,21 @@ def create(
                 index column name is required.
             domain:
                 An optional sequence of tuples specifying the domain of each
-                index column. Each tuple should be a pair consisting of the minimum and
-                maximum values storable in the index column. For example, if there is a
-                single int64-valued index column, then ``domain`` might be ``[(100,
-                200)]`` to indicate that values between 100 and 200, inclusive, can be
-                stored in that column.  If provided, this sequence must have the same
-                length as ``index_column_names``, and the index-column domain will be as
-                specified.  If omitted entirely, or if ``None`` in a given dimension,
-                the corresponding index-column domain will use the minimum and maximum
-                possible values for the column's datatype.  This makes a
-                :class:`DataFrame` growable.
+                index column. Each tuple must be a pair consisting of the
+                minimum and maximum values storable in the index column. For
+                example, if there is a single int64-valued index column, then
+                ``domain`` might be ``[(100, 200)]`` to indicate that values
+                between 100 and 200, inclusive, can be stored in that column.
+                If provided, this sequence must have the same length as
+                ``index_column_names``, and the index-column domain will be as
+                specified.  If omitted entirely, or if ``None`` in a given
+                dimension, the corresponding index-column domain will use an
+                empty range, and data writes after that will fail with "A range
+                was set outside of the current domain". Unless you have a
+                particular reason not to, you should always provide the desired
+                `domain` at create time: this is an optional but strongly
+                recommended parameter. See also ``change_domain`` which allows
+                you to expand the domain after create.
             platform_config:
                 Platform-specific options used to create this array.
                 This may be provided as settings in a dictionary, with options
@@ -621,6 +625,9 @@ def change_domain(
         Lastly, it is an error to try to set the ``domain`` to be smaller than
         ``maxdomain`` along any index column.  The ``maxdomain`` of a dataframe is
         set at creation time, and cannot be extended afterward.
+
+        Lifecycle:
+            Maturing.
         """
         frame = inspect.currentframe()
         function_name_for_messages = frame.f_code.co_name if frame else "tiledbsoma"
@@ -714,31 +721,16 @@ def read(
         _util.check_unpartitioned(partitions)
         self._check_open_read()
 
-        handle = self._handle._handle
-
-        context = handle.context()
-        if platform_config is not None:
-            config = context.tiledb_config.copy()
-            config.update(platform_config)
-            context = clib.SOMAContext(config)
-
-        sr = clib.SOMADataFrame.open(
-            uri=handle.uri,
-            mode=clib.OpenMode.read,
-            context=context,
-            column_names=column_names or [],
+        # TODO: batch_size
+        return TableReadIter(
+            array=self,
+            coords=coords,
+            column_names=column_names,
             result_order=_util.to_clib_result_order(result_order),
-            timestamp=handle.timestamp and (0, handle.timestamp),
+            value_filter=value_filter,
+            platform_config=platform_config,
         )
 
-        if value_filter is not None:
-            sr.set_condition(QueryCondition(value_filter), handle.schema)
-
-        _util._set_coords(sr, coords)
-
-        # TODO: batch_size
-        return TableReadIter(sr)
-
     def write(
         self, values: pa.Table, platform_config: Optional[options.PlatformConfig] = None
     ) -> Self:
diff --git a/apis/python/src/tiledbsoma/_dense_nd_array.py b/apis/python/src/tiledbsoma/_dense_nd_array.py
index a5c9334c7c..ef65f7c9d9 100644
--- a/apis/python/src/tiledbsoma/_dense_nd_array.py
+++ b/apis/python/src/tiledbsoma/_dense_nd_array.py
@@ -20,8 +20,9 @@
 from ._arrow_types import pyarrow_to_carrow_type
 from ._common_nd_array import NDArray
 from ._exception import SOMAError, map_exception_for_create
+from ._read_iters import TableReadIter
 from ._tdb_handles import DenseNDArrayWrapper
-from ._types import OpenTimestamp, Slice, StatusAndReason
+from ._types import OpenTimestamp, Slice
 from ._util import dense_indices_to_shape
 from .options._soma_tiledb_context import (
     SOMATileDBContext,
@@ -232,42 +233,20 @@ def read(
         data_shape = tuple(handle.shape if use_shape else ned)
         target_shape = dense_indices_to_shape(coords, data_shape, result_order)
 
-        context = handle.context()
-        if platform_config is not None:
-            config = context.tiledb_config.copy()
-            config.update(platform_config)
-            context = clib.SOMAContext(config)
-
-        sr = clib.SOMADenseNDArray.open(
-            uri=handle.uri,
-            mode=clib.OpenMode.read,
-            context=context,
+        arrow_table = TableReadIter(
+            array=self,
+            coords=coords,
             column_names=[],
             result_order=_util.to_clib_result_order(result_order),
-            timestamp=handle.timestamp and (0, handle.timestamp),
-        )
-
-        _util._set_coords(sr, coords)
+            value_filter=None,
+            platform_config=platform_config,
+        ).concat()
 
-        arrow_tables = []
-        while True:
-            arrow_table_piece = sr.read_next()
-            if not arrow_table_piece:
-                break
-            arrow_tables.append(arrow_table_piece)
-
-        # For dense arrays there is no zero-output case: attempting to make a test case
-        # to do that, say by indexing a 10x20 array by positions 888 and 999, results
-        # in read-time errors of the form
-        #
-        # [TileDB::Subarray] Error: Cannot add range to dimension 'soma_dim_0'; Range [888, 888] is
-        # out of domain bounds [0, 9]
-        if not arrow_tables:
+        if arrow_table is None:
             raise SOMAError(
                 "internal error: at least one table-piece should have been returned"
             )
 
-        arrow_table = pa.concat_tables(arrow_tables)
         npval = arrow_table.column("soma_data").to_numpy()
         # TODO: as currently coded we're looking at the non-empty domain upper
         # bound but not its lower bound. That works fine if data are written at
@@ -310,7 +289,7 @@ def write(
         """
         _util.check_type("values", values, (pa.Tensor,))
 
-        clib_dense_array = self._handle._handle
+        clib_handle = self._handle._handle
 
         # Compute the coordinates for the dense array.
         new_coords: List[Union[int, Slice[int], None]] = []
@@ -331,34 +310,18 @@ def write(
             if not input.flags.contiguous:
                 input = np.ascontiguousarray(input)
             order = clib.ResultOrder.rowmajor
-        clib_dense_array.reset(result_order=order)
-        _util._set_coords(clib_dense_array, new_coords)
-        clib_dense_array.write(input)
+
+        mq = clib.ManagedQuery(clib_handle, clib_handle.context())
+        mq.set_layout(order)
+        _util._set_coords(mq, clib_handle, new_coords)
+        mq.set_soma_data(input)
+        mq.submit_write()
 
         tiledb_write_options = TileDBWriteOptions.from_platform_config(platform_config)
         if tiledb_write_options.consolidate_and_vacuum:
-            clib_dense_array.consolidate_and_vacuum()
+            clib_handle.consolidate_and_vacuum()
         return self
 
-    def resize(self, newshape: Sequence[Union[int, None]]) -> None:
-        """Supported for ``SparseNDArray``; scheduled for implementation for
-        ``DenseNDArray`` in TileDB-SOMA 1.15
-        """
-        self._handle.resize(newshape)
-
-    def tiledbsoma_upgrade_shape(
-        self, newshape: Sequence[Union[int, None]], check_only: bool = False
-    ) -> StatusAndReason:
-        """Allows the array to have a resizeable shape as described in the TileDB-SOMA
-        1.15 release notes.  Raises an error if the new shape exceeds maxshape in
-        any dimension. Raises an error if the array already has a shape.
-        """
-        if check_only:
-            return self._handle.tiledbsoma_can_upgrade_shape(newshape)
-        else:
-            self._handle.tiledbsoma_upgrade_shape(newshape)
-            return (True, "")
-
     @classmethod
     def _dim_capacity_and_extent(
         cls,
diff --git a/apis/python/src/tiledbsoma/_point_cloud_dataframe.py b/apis/python/src/tiledbsoma/_point_cloud_dataframe.py
index fdfe4fc241..df5bd44baa 100644
--- a/apis/python/src/tiledbsoma/_point_cloud_dataframe.py
+++ b/apis/python/src/tiledbsoma/_point_cloud_dataframe.py
@@ -28,7 +28,6 @@
     _revise_domain_for_extent,
 )
 from ._exception import SOMAError, map_exception_for_create
-from ._query_condition import QueryCondition
 from ._read_iters import TableReadIter
 from ._spatial_dataframe import SpatialDataFrame
 from ._spatial_util import (
@@ -95,13 +94,22 @@ def create(
                 implementation, a ValueError will be raised.
             coordinate_space: Either the coordinate space or the axis names for the
                 coordinate space the point cloud is defined on.
-            domain: An optional sequence of tuples specifying the domain of each
-                index column. Each tuple should be a pair consisting of the minimum
-                and maximum values storable in the index column. If omitted entirely,
-                or if ``None`` in a given dimension, the corresponding index-column
-                domain will use the minimum and maximum possible values for the
-                column's datatype.  This makes a point cloud dataframe growable.
-
+            domain:
+                An optional sequence of tuples specifying the domain of each
+                index column. Each tuple must be a pair consisting of the
+                minimum and maximum values storable in the index column. For
+                example, if there is a single int64-valued index column, then
+                ``domain`` might be ``[(100, 200)]`` to indicate that values
+                between 100 and 200, inclusive, can be stored in that column.
+                If provided, this sequence must have the same length as
+                ``index_column_names``, and the index-column domain will be as
+                specified.  If omitted entirely, or if ``None`` in a given
+                dimension, the corresponding index-column domain will use an
+                empty range, and data writes after that will fail with "A range
+                was set outside of the current domain". Unless you have a
+                particular reason not to, you should always provide the desired
+                `domain` at create time: this is an optional but strongly
+                recommended parameter.
         Returns:
             The newly created point cloud, opened for writing.
 
@@ -323,31 +331,16 @@ def read(
         _util.check_unpartitioned(partitions)
         self._check_open_read()
 
-        handle = self._handle._handle
-
-        context = handle.context()
-        if platform_config is not None:
-            config = context.tiledb_config.copy()
-            config.update(platform_config)
-            context = clib.SOMAContext(config)
-
-        sr = clib.SOMAPointCloudDataFrame.open(
-            uri=handle.uri,
-            mode=clib.OpenMode.read,
-            context=context,
-            column_names=column_names or [],
+        # TODO: batch_size
+        return TableReadIter(
+            array=self,
+            coords=coords,
+            column_names=column_names,
             result_order=_util.to_clib_result_order(result_order),
-            timestamp=handle.timestamp and (0, handle.timestamp),
+            value_filter=value_filter,
+            platform_config=platform_config,
         )
 
-        if value_filter is not None:
-            sr.set_condition(QueryCondition(value_filter), handle.schema)
-
-        _util._set_coords(sr, coords)
-
-        # # TODO: batch_size
-        return TableReadIter(sr)
-
     def read_spatial_region(
         self,
         region: Optional[options.SpatialRegion] = None,
diff --git a/apis/python/src/tiledbsoma/_read_iters.py b/apis/python/src/tiledbsoma/_read_iters.py
index bec8cc81d8..92a8198ca7 100644
--- a/apis/python/src/tiledbsoma/_read_iters.py
+++ b/apis/python/src/tiledbsoma/_read_iters.py
@@ -13,7 +13,6 @@
 from typing import (
     TYPE_CHECKING,
     Any,
-    Dict,
     Iterator,
     List,
     Optional,
@@ -39,11 +38,13 @@
 from ._exception import SOMAError
 from ._fastercsx import CompressedMatrix
 from ._indexer import IntIndexer
+from ._query_condition import QueryCondition
 from ._types import NTuple
 from .options import SOMATileDBContext
 
 if TYPE_CHECKING:
     from . import SparseNDArray
+    from ._soma_array import SOMAArray
 
 
 # Convenience types
@@ -66,8 +67,47 @@
 class TableReadIter(somacore.ReadIter[pa.Table]):
     """Iterator over `Arrow Table <https://arrow.apache.org/docs/python/generated/pyarrow.Table.html>`_ elements"""
 
-    def __init__(self, sr: clib.SOMAArray):
-        self._reader = _arrow_table_reader(sr)
+    def __init__(
+        self,
+        array: SOMAArray,
+        coords: Union[
+            options.SparseDFCoords, options.SparseNDCoords, options.DenseNDCoords
+        ],
+        column_names: Optional[Sequence[str]],
+        result_order: clib.ResultOrder,
+        value_filter: Optional[str],
+        platform_config: Optional[options.PlatformConfig],
+    ):
+        """Initalizes a new TableReadIter for SOMAArrays.
+
+        Args:
+            array (SOMAArray):
+                The NDArray, DataFrame, or SpatialDataFrame being read.
+
+            coords (Union[
+                options.SparseDFCoords, options.SparseNDCoords, options.DenseNDCoords
+            ]):
+                for each index dimension, which rows to read.
+                ``()`` means no constraint -- all IDs.
+
+            column_names (Optional[Sequence[str]]):
+                The named columns to read and return.
+                ``None`` means no constraint -- all column names.
+
+            result_order (clib.ResultOrder):
+                Order of read results.
+                This can be one of automatic, rowmajor, or colmajor.
+
+            value_filter (Optional[str]):
+                An optional [value filter] to apply to the results.
+
+            platform_config (Optional[options.PlatformConfig]):
+                Pass in parameters for tuning reads.
+
+        """
+        self._reader = ArrowTableRead(
+            array, coords, column_names, result_order, value_filter, platform_config
+        )
 
     def __next__(self) -> pa.Table:
         return next(self._reader)
@@ -87,10 +127,11 @@ class BlockwiseReadIterBase(somacore.ReadIter[_RT], metaclass=abc.ABCMeta):
 
     def __init__(
         self,
-        array: "SparseNDArray",
-        sr: clib.SOMAArray,
+        array: SOMAArray,
         coords: options.SparseNDCoords,
         axis: Union[int, Sequence[int]],
+        result_order: clib.ResultOrder,
+        platform_config: Optional[options.PlatformConfig],
         *,
         size: Optional[Union[int, Sequence[int]]] = None,
         reindex_disable_on_axis: Optional[Union[int, Sequence[int]]] = None,
@@ -99,10 +140,12 @@ def __init__(
     ):
         super().__init__()
 
-        self.ndim = len(sr.shape)
         self.array = array
-        self.sr = sr
+        self.shape = array._handle._handle.shape
+        self.ndim = len(self.shape)
         self.eager = eager
+        self.result_order = result_order
+        self.platform_config = platform_config
 
         # Assign a thread pool from the context, or create a new one if no context
         # is available
@@ -115,7 +158,7 @@ def __init__(
 
         # raises on various error checks, AND normalizes args
         self.axis, self.size, self.reindex_disable_on_axis = self._validate_args(
-            sr.shape, axis, size, reindex_disable_on_axis
+            self.shape, axis, size, reindex_disable_on_axis
         )
 
         self.major_axis = self.axis[0]
@@ -125,11 +168,7 @@ def __init__(
         self.joinids: List[pa.Array] = [
             pa.array(
                 np.concatenate(
-                    list(
-                        _coords_strider(
-                            self.coords[d], self.sr.shape[d], self.sr.shape[d]
-                        )
-                    )
+                    list(_coords_strider(self.coords[d], self.shape[d], self.shape[d]))
                 )
                 if d != self.major_axis
                 else np.array([], dtype=np.int64)
@@ -230,20 +269,26 @@ def _maybe_eager_iterator(
 
     def _table_reader(self) -> Iterator[BlockwiseTableReadIterResult]:
         """Private. Blockwise table reader. Helper function for sub-class use"""
-        kwargs: Dict[str, object] = {"result_order": self.sr.result_order}
         for coord_chunk in _coords_strider(
             self.coords[self.major_axis],
-            self.sr.shape[self.major_axis],
+            self.shape[self.major_axis],
             self.size[0],
         ):
-            self.sr.reset(**kwargs)
             step_coords = list(self.coords)
             step_coords[self.major_axis] = coord_chunk
-            _util._set_coords(self.sr, step_coords)
 
             joinids = list(self.joinids)
             joinids[self.major_axis] = pa.array(coord_chunk)
-            yield pa.concat_tables(_arrow_table_reader(self.sr)), tuple(joinids)
+            yield pa.concat_tables(
+                ArrowTableRead(
+                    array=self.array,
+                    coords=step_coords,
+                    column_names=[],  # select all columns
+                    result_order=self.result_order,
+                    value_filter=None,
+                    platform_config=self.platform_config,
+                )
+            ), tuple(joinids)
 
     def _reindexed_table_reader(
         self,
@@ -286,10 +331,11 @@ class BlockwiseScipyReadIter(BlockwiseReadIterBase[BlockwiseScipyReadIterResult]
 
     def __init__(
         self,
-        array: "SparseNDArray",
-        sr: clib.SOMAArray,
+        array: SOMAArray,
         coords: options.SparseNDCoords,
         axis: Union[int, Sequence[int]],
+        result_order: clib.ResultOrder,
+        platform_config: Optional[options.PlatformConfig],
         *,
         size: Optional[Union[int, Sequence[int]]] = None,
         reindex_disable_on_axis: Optional[Union[int, Sequence[int]]] = None,
@@ -301,9 +347,10 @@ def __init__(
         self.context = context
         super().__init__(
             array,
-            sr,
             coords,
             axis,
+            result_order,
+            platform_config,
             size=size,
             reindex_disable_on_axis=reindex_disable_on_axis,
             eager=eager,
@@ -311,7 +358,7 @@ def __init__(
         )
 
         if (
-            len(self.sr.shape) != 2
+            len(self.shape) != 2
             or len(self.coords) > 2
             or self.major_axis not in [0, 1]
         ):
@@ -365,7 +412,7 @@ def _mk_shape(
         self, major_coords: npt.NDArray[np.int64], minor_coords: npt.NDArray[np.int64]
     ) -> Tuple[int, int]:
         """Private. Make shape of this iterator step"""
-        shape = cast(Tuple[int, int], tuple(self.sr.shape))
+        shape = cast(Tuple[int, int], tuple(self.shape))
         assert len(shape) == 2
         _sp_shape: List[int] = list(shape)
 
@@ -393,7 +440,7 @@ def _coo_reader(
             )
 
             # SOMA disallows duplicates. Canonical implies sorted row-major, no dups
-            if self.sr.result_order == clib.ResultOrder.rowmajor:
+            if self.result_order == clib.ResultOrder.rowmajor:
                 sp.has_canonical_format = True
 
             yield sp, indices
@@ -426,20 +473,41 @@ def _cs_reader(
 class SparseTensorReadIterBase(somacore.ReadIter[_RT], metaclass=abc.ABCMeta):
     """Private implementation class"""
 
-    def __init__(self, sr: clib.SOMAArray, shape: NTuple):
-        self.sr = sr
+    def __init__(
+        self,
+        array: SparseNDArray,
+        coords: options.SparseDFCoords,
+        shape: NTuple,
+        result_order: clib.ResultOrder,
+        platform_config: Optional[options.PlatformConfig],
+    ):
+        self.array = array
+        self.coords = coords
         self.shape = shape
+        self.result_order = result_order
+        self.platform_config = platform_config
+
+        clib_handle = array._handle._handle
+
+        if platform_config is not None:
+            cfg = clib_handle.context().config()
+            cfg.update(platform_config)
+            ctx = clib.SOMAContext(cfg)
+        else:
+            ctx = clib_handle.context()
+
+        self.mq = clib.ManagedQuery(clib_handle, ctx)
+
+        self.mq.set_layout(result_order)
+
+        _util._set_coords(self.mq, clib_handle, coords)
 
     @abc.abstractmethod
     def _from_table(self, arrow_table: pa.Table) -> _RT:
         raise NotImplementedError()
 
     def __next__(self) -> _RT:
-        arrow_table = self.sr.read_next()
-        if arrow_table is None:
-            raise StopIteration
-
-        return self._from_table(arrow_table)
+        return self._from_table(self.mq.next())
 
     def concat(self) -> _RT:
         """Returns all the requested data in a single operation.
@@ -447,7 +515,16 @@ def concat(self) -> _RT:
         If some data has already been retrieved using ``next``, this will return
         the rest of the data after that is already returned.
         """
-        arrow_tables = pa.concat_tables(TableReadIter(self.sr))
+        arrow_tables = pa.concat_tables(
+            TableReadIter(
+                array=self.array,
+                coords=self.coords,
+                column_names=[],  # select all columns
+                result_order=self.result_order,
+                value_filter=None,
+                platform_config=self.platform_config,
+            )
+        )
         return self._from_table(arrow_tables)
 
 
@@ -465,12 +542,41 @@ def _from_table(self, arrow_table: pa.Table) -> pa.SparseCOOTensor:
         return pa.SparseCOOTensor.from_numpy(coo_data, coo_coords, shape=self.shape)
 
 
-def _arrow_table_reader(sr: clib.SOMAArray) -> Iterator[pa.Table]:
-    """Private. Simple Table iterator on any Array"""
-    tbl = sr.read_next()
-    while tbl is not None:
-        yield tbl
-        tbl = sr.read_next()
+class ArrowTableRead(Iterator[pa.Table]):
+    def __init__(
+        self,
+        array: SOMAArray,
+        coords: Union[
+            options.SparseDFCoords, options.SparseNDCoords, options.DenseNDCoords
+        ],
+        column_names: Optional[Sequence[str]],
+        result_order: clib.ResultOrder,
+        value_filter: Optional[str],
+        platform_config: Optional[options.PlatformConfig],
+    ):
+        clib_handle = array._handle._handle
+
+        if platform_config is not None:
+            cfg = clib_handle.context().config()
+            cfg.update(platform_config)
+            ctx = clib.SOMAContext(cfg)
+        else:
+            ctx = clib_handle.context()
+
+        self.mq = clib.ManagedQuery(clib_handle, ctx)
+
+        self.mq.set_layout(result_order)
+
+        if column_names is not None:
+            self.mq.select_columns(list(column_names))
+
+        if value_filter is not None:
+            self.mq.set_condition(QueryCondition(value_filter), clib_handle.schema)
+
+        _util._set_coords(self.mq, clib_handle, coords)
+
+    def __next__(self) -> pa.Table:
+        return self.mq.next()
 
 
 def _coords_strider(
diff --git a/apis/python/src/tiledbsoma/_scene.py b/apis/python/src/tiledbsoma/_scene.py
index f7432c3f73..71aabdda27 100644
--- a/apis/python/src/tiledbsoma/_scene.py
+++ b/apis/python/src/tiledbsoma/_scene.py
@@ -6,6 +6,7 @@
 Implementation of a SOMA Scene
 """
 
+import warnings
 from typing import Any, List, Optional, Sequence, Tuple, Type, TypeVar, Union
 
 import somacore
@@ -17,9 +18,13 @@
 from typing_extensions import Self
 
 from . import _funcs, _tdb_handles
+from . import pytiledbsoma as clib
 from ._collection import CollectionBase
-from ._constants import SOMA_COORDINATE_SPACE_METADATA_KEY
-from ._exception import SOMAError
+from ._constants import (
+    SOMA_COORDINATE_SPACE_METADATA_KEY,
+    SPATIAL_DISCLAIMER,
+)
+from ._exception import SOMAError, map_exception_for_create
 from ._geometry_dataframe import GeometryDataFrame
 from ._multiscale_image import MultiscaleImage
 from ._point_cloud_dataframe import PointCloudDataFrame
@@ -32,6 +37,7 @@
 )
 from ._types import OpenTimestamp
 from .options import SOMATileDBContext
+from .options._soma_tiledb_context import _validate_soma_tiledb_context
 
 _spatial_element = Union[GeometryDataFrame, MultiscaleImage, PointCloudDataFrame]
 
@@ -97,16 +103,29 @@ def create(
         Lifecycle:
             Experimental.
         """
-        if coordinate_space is not None:
-            raise NotImplementedError(
-                "Setting the coordinate space on create is not yet implemented."
+        warnings.warn(SPATIAL_DISCLAIMER)
+
+        context = _validate_soma_tiledb_context(context)
+        try:
+            timestamp_ms = context._open_timestamp_ms(tiledb_timestamp)
+            clib.SOMAScene.create(
+                ctx=context.native_context,
+                uri=uri,
+                timestamp=(0, timestamp_ms),
             )
-        return super().create(
-            uri=uri,
-            platform_config=platform_config,
-            context=context,
-            tiledb_timestamp=tiledb_timestamp,
-        )
+            handle = cls._wrapper_type.open(uri, "w", context, tiledb_timestamp)
+            if coordinate_space is not None:
+                if not isinstance(coordinate_space, CoordinateSpace):
+                    coordinate_space = CoordinateSpace.from_axis_names(coordinate_space)
+                handle.meta[SOMA_COORDINATE_SPACE_METADATA_KEY] = (
+                    coordinate_space_to_json(coordinate_space)
+                )
+            return cls(
+                handle,
+                _dont_call_this_use_create_or_open_instead="tiledbsoma-internal-code",
+            )
+        except SOMAError as e:
+            raise map_exception_for_create(e, uri) from None
 
     def __init__(
         self,
diff --git a/apis/python/src/tiledbsoma/_sparse_nd_array.py b/apis/python/src/tiledbsoma/_sparse_nd_array.py
index 55ab1678f1..60891c9629 100644
--- a/apis/python/src/tiledbsoma/_sparse_nd_array.py
+++ b/apis/python/src/tiledbsoma/_sparse_nd_array.py
@@ -40,7 +40,7 @@
     TableReadIter,
 )
 from ._tdb_handles import SparseNDArrayWrapper
-from ._types import NTuple, OpenTimestamp, StatusAndReason
+from ._types import NTuple, OpenTimestamp
 from .options._soma_tiledb_context import (
     SOMATileDBContext,
     _validate_soma_tiledb_context,
@@ -266,58 +266,16 @@ def read(
             * Negative indexing is unsupported.
         """
         del batch_size  # Currently unused.
-        handle: clib.SOMASparseNDArray = self._handle._handle
-
         self._check_open_read()
         _util.check_unpartitioned(partitions)
 
-        context = handle.context()
-        if platform_config is not None:
-            config = context.tiledb_config.copy()
-            config.update(platform_config)
-            context = clib.SOMAContext(config)
-
-        sr = clib.SOMASparseNDArray.open(
-            uri=handle.uri,
-            mode=clib.OpenMode.read,
-            context=context,
-            column_names=[],
+        return SparseNDArrayRead(
+            array=self,
+            coords=coords,
             result_order=_util.to_clib_result_order(result_order),
-            timestamp=handle.timestamp and (0, handle.timestamp),
+            platform_config=platform_config,
         )
 
-        return SparseNDArrayRead(sr, self, coords)
-
-    def resize(
-        self, newshape: Sequence[Union[int, None]], check_only: bool = False
-    ) -> StatusAndReason:
-        """Increases the shape of the array as specfied. Raises an error if the new
-        shape is less than the current shape in any dimension. Raises an error if
-        the new shape exceeds maxshape in any dimension. Raises an error if the
-        array doesn't already have a shape: in that case please call
-        tiledbsoma_upgrade_shape. If ``check_only`` is ``True``, returns
-        whether the operation would succeed if attempted, and a reason why it
-        would not.
-        """
-        if check_only:
-            return self._handle.tiledbsoma_can_resize(newshape)
-        else:
-            self._handle.resize(newshape)
-            return (True, "")
-
-    def tiledbsoma_upgrade_shape(
-        self, newshape: Sequence[Union[int, None]], check_only: bool = False
-    ) -> StatusAndReason:
-        """Allows the array to have a resizeable shape as described in the TileDB-SOMA
-        1.15 release notes.  Raises an error if the new shape exceeds maxshape in
-        any dimension. Raises an error if the array already has a shape.
-        """
-        if check_only:
-            return self._handle.tiledbsoma_can_upgrade_shape(newshape)
-        else:
-            self._handle.tiledbsoma_upgrade_shape(newshape)
-            return (True, "")
-
     def write(
         self,
         values: Union[
@@ -555,18 +513,20 @@ class _SparseNDArrayReadBase(somacore.SparseRead):
 
     def __init__(
         self,
-        sr: clib.SOMAArray,
         array: SparseNDArray,
         coords: options.SparseNDCoords,
+        result_order: clib.ResultOrder,
+        platform_config: Optional[options.PlatformConfig],
     ):
         """
         Lifecycle:
             Maturing.
         """
-        self.sr = sr
-        self.shape = tuple(sr.shape)
         self.array = array
         self.coords = coords
+        self.shape = tuple(array._handle._handle.shape)
+        self.result_order = result_order
+        self.platform_config = platform_config
 
 
 class SparseNDArrayRead(_SparseNDArrayReadBase):
@@ -600,8 +560,13 @@ def coos(self, shape: Optional[NTuple] = None) -> SparseCOOTensorReadIter:
         """
         if shape is not None and (len(shape) != len(self.shape)):
             raise ValueError(f"shape must be a tuple of size {len(self.shape)}")
-        _util._set_coords(self.sr, self.coords)
-        return SparseCOOTensorReadIter(self.sr, shape or self.shape)
+        return SparseCOOTensorReadIter(
+            self.array,
+            self.coords,
+            shape or self.shape,
+            self.result_order,
+            self.platform_config,
+        )
 
     def tables(self) -> TableReadIter:
         """
@@ -611,8 +576,14 @@ def tables(self) -> TableReadIter:
         Lifecycle:
             Maturing.
         """
-        _util._set_coords(self.sr, self.coords)
-        return TableReadIter(self.sr)
+        return TableReadIter(
+            array=self.array,
+            coords=self.coords,
+            column_names=[],
+            result_order=self.result_order,
+            value_filter=None,
+            platform_config=self.platform_config,
+        )
 
     def blockwise(
         self,
@@ -683,10 +654,11 @@ def blockwise(
             Maturing.
         """
         return SparseNDArrayBlockwiseRead(
-            self.sr,
             self.array,
             self.coords,
             axis,
+            self.result_order,
+            self.platform_config,
             size=size,
             reindex_disable_on_axis=reindex_disable_on_axis,
             eager=eager,
@@ -696,16 +668,18 @@ def blockwise(
 class SparseNDArrayBlockwiseRead(_SparseNDArrayReadBase):
     def __init__(
         self,
-        sr: clib.SOMAArray,
         array: SparseNDArray,
         coords: options.SparseNDCoords,
         axis: Union[int, Sequence[int]],
+        result_order: clib.ResultOrder,
+        platform_config: Optional[options.PlatformConfig],
         *,
         size: Optional[Union[int, Sequence[int]]],
         reindex_disable_on_axis: Optional[Union[int, Sequence[int]]],
         eager: bool = True,
     ):
-        super().__init__(sr, array, coords)
+        super().__init__(array, coords, result_order, platform_config)
+        self.result_order = result_order
         self.axis = axis
         self.size = size
         self.reindex_disable_on_axis = reindex_disable_on_axis
@@ -726,9 +700,10 @@ def tables(self) -> BlockwiseTableReadIter:
         """
         return BlockwiseTableReadIter(
             self.array,
-            self.sr,
             self.coords,
             self.axis,
+            self.result_order,
+            self.platform_config,
             size=self.size,
             reindex_disable_on_axis=self.reindex_disable_on_axis,
             eager=self.eager,
@@ -771,9 +746,10 @@ def scipy(self, *, compress: bool = True) -> BlockwiseScipyReadIter:
         """
         return BlockwiseScipyReadIter(
             self.array,
-            self.sr,
             self.coords,
             self.axis,
+            self.result_order,
+            self.platform_config,
             size=self.size,
             compress=compress,
             reindex_disable_on_axis=self.reindex_disable_on_axis,
diff --git a/apis/python/src/tiledbsoma/_util.py b/apis/python/src/tiledbsoma/_util.py
index ea118f82e0..e007ab9c16 100644
--- a/apis/python/src/tiledbsoma/_util.py
+++ b/apis/python/src/tiledbsoma/_util.py
@@ -464,7 +464,9 @@ def _cast_domainish(domainish: List[Any]) -> Tuple[Tuple[object, object], ...]:
     return tuple(result)
 
 
-def _set_coords(sarr: clib.SOMAArray, coords: options.SparseNDCoords) -> None:
+def _set_coords(
+    mq: clib.ManagedQuery, sarr: clib.SOMAArray, coords: options.SparseNDCoords
+) -> None:
     if not is_nonstringy_sequence(coords):
         raise TypeError(
             f"coords type {type(coords)} must be a regular sequence,"
@@ -478,10 +480,12 @@ def _set_coords(sarr: clib.SOMAArray, coords: options.SparseNDCoords) -> None:
         )
 
     for i, coord in enumerate(coords):
-        _set_coord(i, sarr, coord)
+        _set_coord(i, mq, sarr, coord)
 
 
-def _set_coord(dim_idx: int, sarr: clib.SOMAArray, coord: object) -> None:
+def _set_coord(
+    dim_idx: int, mq: clib.ManagedQuery, sarr: clib.SOMAArray, coord: object
+) -> None:
     if coord is None:
         return
 
@@ -489,19 +493,19 @@ def _set_coord(dim_idx: int, sarr: clib.SOMAArray, coord: object) -> None:
     dom = _cast_domainish(sarr.domain())[dim_idx]
 
     if isinstance(coord, (str, bytes)):
-        sarr.set_dim_points_string_or_bytes(dim.name, [coord])
+        mq.set_dim_points_string_or_bytes(dim.name, [coord])
         return
 
     if isinstance(coord, (pa.Array, pa.ChunkedArray)):
-        sarr.set_dim_points_arrow(dim.name, coord)
+        mq.set_dim_points_arrow(dim.name, coord)
         return
 
     if isinstance(coord, (Sequence, np.ndarray)):
-        _set_coord_by_py_seq_or_np_array(sarr, dim, coord)
+        _set_coord_by_py_seq_or_np_array(mq, dim, coord)
         return
 
     if isinstance(coord, int):
-        sarr.set_dim_points_int64(dim.name, [coord])
+        mq.set_dim_points_int64(dim.name, [coord])
         return
 
     # Note: slice(None, None) matches the is_slice_of part, unless we also check
@@ -520,7 +524,7 @@ def _set_coord(dim_idx: int, sarr: clib.SOMAArray, coord: object) -> None:
             _, stop = ned[dim_idx]
         else:
             stop = coord.stop
-        sarr.set_dim_ranges_string_or_bytes(dim.name, [(start, stop)])
+        mq.set_dim_ranges_string_or_bytes(dim.name, [(start, stop)])
         return
 
     # Note: slice(None, None) matches the is_slice_of part, unless we also check
@@ -543,21 +547,21 @@ def _set_coord(dim_idx: int, sarr: clib.SOMAArray, coord: object) -> None:
         else:
             istop = ts_dom[1].as_py()
 
-        sarr.set_dim_ranges_int64(dim.name, [(istart, istop)])
+        mq.set_dim_ranges_int64(dim.name, [(istart, istop)])
         return
 
     if isinstance(coord, slice):
         validate_slice(coord)
         if coord.start is None and coord.stop is None:
             return
-        _set_coord_by_numeric_slice(sarr, dim, dom, coord)
+        _set_coord_by_numeric_slice(mq, dim, dom, coord)
         return
 
     raise TypeError(f"unhandled type {dim.type} for index column named {dim.name}")
 
 
 def _set_coord_by_py_seq_or_np_array(
-    sarr: clib.SOMAArray, dim: pa.Field, coord: object
+    mq: clib.ManagedQuery, dim: pa.Field, coord: object
 ) -> None:
     if isinstance(coord, np.ndarray):
         if coord.ndim != 1:
@@ -566,7 +570,7 @@ def _set_coord_by_py_seq_or_np_array(
             )
 
     try:
-        set_dim_points = getattr(sarr, f"set_dim_points_{dim.type}")
+        set_dim_points = getattr(mq, f"set_dim_points_{dim.type}")
     except AttributeError:
         # We have to handle this type specially below
         pass
@@ -575,7 +579,7 @@ def _set_coord_by_py_seq_or_np_array(
         return
 
     if pa_types_is_string_or_bytes(dim.type):
-        sarr.set_dim_points_string_or_bytes(dim.name, coord)
+        mq.set_dim_points_string_or_bytes(dim.name, coord)
         return
 
     if pa.types.is_timestamp(dim.type):
@@ -586,14 +590,14 @@ def _set_coord_by_py_seq_or_np_array(
         icoord = [
             int(e.astype("int64")) if isinstance(e, np.datetime64) else e for e in coord
         ]
-        sarr.set_dim_points_int64(dim.name, icoord)
+        mq.set_dim_points_int64(dim.name, icoord)
         return
 
     raise ValueError(f"unhandled type {dim.type} for index column named {dim.name}")
 
 
 def _set_coord_by_numeric_slice(
-    sarr: clib.SOMAArray, dim: pa.Field, dom: Tuple[object, object], coord: Slice[Any]
+    mq: clib.ManagedQuery, dim: pa.Field, dom: Tuple[object, object], coord: Slice[Any]
 ) -> None:
     try:
         lo_hi = slice_to_numeric_range(coord, dom)
@@ -604,7 +608,7 @@ def _set_coord_by_numeric_slice(
         return
 
     try:
-        set_dim_range = getattr(sarr, f"set_dim_ranges_{dim.type}")
+        set_dim_range = getattr(mq, f"set_dim_ranges_{dim.type}")
         set_dim_range(dim.name, [lo_hi])
         return
     except AttributeError:
diff --git a/apis/python/src/tiledbsoma/io/conversions.py b/apis/python/src/tiledbsoma/io/conversions.py
index 107a83855b..0f7380f3a4 100644
--- a/apis/python/src/tiledbsoma/io/conversions.py
+++ b/apis/python/src/tiledbsoma/io/conversions.py
@@ -13,6 +13,7 @@
 import numpy as np
 import pandas as pd
 import pandas._typing as pdt
+import pandas.api.types
 import pyarrow as pa
 import scipy.sparse as sp
 
@@ -25,19 +26,31 @@
 _MT = TypeVar("_MT", NPNDArray, sp.spmatrix, PDSeries)
 _str_to_type = {"boolean": bool, "string": str, "bytes": bytes}
 
+STRING_DECAT_THRESHOLD = 4096
+"""
+For enum-of-string columns with a cardinality higher than this, we convert from
+enum-of-string in the AnnData ``obs``/``var``, to plain string in TileDB-SOMA
+``obs``/``var``. However, if we're appending to existing storage, we follow the
+schema there.
+"""
 
-def decategoricalize_obs_or_var(obs_or_var: pd.DataFrame) -> pd.DataFrame:
-    """Performs a typecast into types that TileDB can persist."""
-    if len(obs_or_var.columns) > 0:
-        return pd.DataFrame.from_dict(
-            {
-                str(k): to_tiledb_supported_array_type(str(k), v)
-                for k, v in obs_or_var.items()
-            },
-        )
-    else:
+
+def obs_or_var_to_tiledb_supported_array_type(obs_or_var: pd.DataFrame) -> pd.DataFrame:
+    """
+    Performs a typecast into types that TileDB can persist.  This includes, as a
+    performance improvement, converting high-cardinality categorical-of-string
+    columns (cardinality > 4096) to plain string.
+    """
+    if len(obs_or_var.columns) == 0:
         return obs_or_var.copy()
 
+    return pd.DataFrame.from_dict(
+        {
+            str(k): to_tiledb_supported_array_type(str(k), v)
+            for k, v in obs_or_var.items()
+        },
+    )
+
 
 @typeguard_ignore
 def _to_tiledb_supported_dtype(dtype: _DT) -> _DT:
@@ -47,8 +60,8 @@ def _to_tiledb_supported_dtype(dtype: _DT) -> _DT:
 
 
 def to_tiledb_supported_array_type(name: str, x: _MT) -> _MT:
-    """Converts datatypes unrepresentable by TileDB into datatypes it can represent.
-    E.g., float16 -> float32
+    """Converts datatypes unrepresentable by TileDB into datatypes it can represent,
+    e.g., float16 -> float32.
     """
     if isinstance(x, (np.ndarray, sp.spmatrix)) or not isinstance(
         x.dtype, pd.CategoricalDtype
@@ -56,27 +69,17 @@ def to_tiledb_supported_array_type(name: str, x: _MT) -> _MT:
         target_dtype = _to_tiledb_supported_dtype(x.dtype)
         return x if target_dtype == x.dtype else x.astype(target_dtype)
 
-    # categories = x.cat.categories
-    # cat_dtype = categories.dtype
-    # if cat_dtype.kind in ("f", "u", "i"):
-    #     if x.hasnans and cat_dtype.kind == "i":
-    #         raise ValueError(
-    #             f"Categorical column {name!r} contains NaN -- unable to convert to TileDB array."
-    #         )
-    #     # More mysterious spurious mypy errors.
-    #     target_dtype = _to_tiledb_supported_dtype(cat_dtype)  # type: ignore[arg-type]
-    # else:
-    #     # Into the weirdness. See if Pandas can help with edge cases.
-    #     inferred = infer_dtype(categories)
-    #     if x.hasnans and inferred in ("boolean", "bytes"):
-    #         raise ValueError(
-    #             "Categorical array contains NaN -- unable to convert to TileDB array."
-    #         )
-    #     target_dtype = np.dtype(  # type: ignore[assignment]
-    #         _str_to_type.get(inferred, object)
-    #     )
-
-    # return x.astype(target_dtype)
+    # If the column is categorical-of-string of high cardinality, we declare
+    # this is likely a mistake, and it will definitely lead to performance
+    # issues in subsequent processing.
+    if isinstance(x, pd.Series) and isinstance(x.dtype, pd.CategoricalDtype):
+        # Heuristic number
+        if (
+            pandas.api.types.is_string_dtype(x)
+            and len(x.cat.categories) > STRING_DECAT_THRESHOLD
+        ):
+            return x.astype(str)
+
     return x
 
 
diff --git a/apis/python/src/tiledbsoma/io/ingest.py b/apis/python/src/tiledbsoma/io/ingest.py
index d776f7a941..cbf239f157 100644
--- a/apis/python/src/tiledbsoma/io/ingest.py
+++ b/apis/python/src/tiledbsoma/io/ingest.py
@@ -516,7 +516,7 @@ def from_anndata(
     df_uri = _util.uri_joinpath(experiment_uri, "obs")
     with _write_dataframe(
         df_uri,
-        conversions.decategoricalize_obs_or_var(anndata.obs),
+        conversions.obs_or_var_to_tiledb_supported_array_type(anndata.obs),
         id_column_name=obs_id_name,
         axis_mapping=jidmaps.obs_axis,
         **ingest_platform_ctx,
@@ -566,7 +566,7 @@ def from_anndata(
             # MS/meas/VAR
             with _write_dataframe(
                 _util.uri_joinpath(measurement_uri, "var"),
-                conversions.decategoricalize_obs_or_var(anndata.var),
+                conversions.obs_or_var_to_tiledb_supported_array_type(anndata.var),
                 id_column_name=var_id_name,
                 # Layer existence is pre-checked in the registration phase
                 axis_mapping=jidmaps.var_axes[measurement_name],
@@ -701,7 +701,9 @@ def _ingest_obs_var_m_p(
 
                         with _write_dataframe(
                             _util.uri_joinpath(raw_uri, "var"),
-                            conversions.decategoricalize_obs_or_var(anndata.raw.var),
+                            conversions.obs_or_var_to_tiledb_supported_array_type(
+                                anndata.raw.var
+                            ),
                             id_column_name=var_id_name,
                             axis_mapping=jidmaps.var_axes["raw"],
                             **ingest_platform_ctx,
@@ -794,7 +796,7 @@ def append_obs(
 
     with _write_dataframe(
         exp.obs.uri,
-        conversions.decategoricalize_obs_or_var(new_obs),
+        conversions.obs_or_var_to_tiledb_supported_array_type(new_obs),
         id_column_name=obs_id_name,
         platform_config=platform_config,
         context=context,
@@ -859,7 +861,7 @@ def append_var(
 
     with _write_dataframe(
         sdf.uri,
-        conversions.decategoricalize_obs_or_var(new_var),
+        conversions.obs_or_var_to_tiledb_supported_array_type(new_var),
         id_column_name=var_id_name,
         platform_config=platform_config,
         context=context,
@@ -1068,6 +1070,130 @@ def _create_or_open_coll(
     )
 
 
+def _extract_new_values_for_append_aux(
+    previous_soma_dataframe: DataFrame,
+    arrow_table: pa.Table,
+) -> pa.Table:
+    """
+    Helper function for _extract_new_values_for_append.
+
+    This does two things:
+
+    * Retains only the 'new' rows compared to existing storage.
+      Example is append-mode updates to the var dataframe for
+      which it's likely that most/all gene IDs have already been seen.
+
+    * String vs categorical of string:
+
+      o If we're appending a plain-string column to existing
+        categorical-of-string storage, convert the about-to-be-written data
+        to categorical of string, to match.
+
+      o If we're appending a categorical-of-string column to existing
+        plain-string storage, convert the about-to-be-written data
+        to plain string, to match.
+
+    Context: https://github.com/single-cell-data/TileDB-SOMA/issues/3353.
+    Namely, we find that AnnData's to_h5ad/from_h5ad can categoricalize (without
+    the user's knowledge or intention) string columns.  For example, even
+    cell_id/barcode, for which there may be millions of distinct values, with no
+    gain to be had from dictionary encoding, will be converted to categorical.
+    We find that converting these high-cardinality enums to plain string is a
+    significant performance win for subsequent accesses. When we do an initial
+    ingest from AnnData to TileDB-SOMA, we convert from categorical-of-string to
+    plain string if the cardinality exceeds some threshold.
+
+    All well and good -- except for one more complication which is append mode.
+    Namely, if the new column has high enough cardinality that we would
+    downgrade to plain string, but the existing storage has
+    categorical-of-string, we must write the new data as categorical-of-string.
+    Likewise, if the new column has low enough cardinality that we would keep it
+    as categorical-of-string, but the existing storage has plain string, we must
+    write the new data as plain strings.
+    """
+
+    # Retain only the new rows.
+    previous_sjids_table = previous_soma_dataframe.read(
+        column_names=["soma_joinid"]
+    ).concat()
+    previous_join_ids = set(
+        int(e)
+        for e in get_dataframe_values(previous_sjids_table.to_pandas(), SOMA_JOINID)
+    )
+    mask = [e.as_py() not in previous_join_ids for e in arrow_table[SOMA_JOINID]]
+    arrow_table = arrow_table.filter(mask)
+
+    # This is a redundant, failsafe check. The append-mode registrar already
+    # ensure schema homogeneity before we get here.
+    old_schema = previous_soma_dataframe.schema
+    new_schema = arrow_table.schema
+
+    # Note: we may be doing an add-column that doesn't exist in tiledbsoma
+    # storage but is present in the new AnnData. We don't need to change
+    # anything in that case. Regardless, we can't assume that the old
+    # and new schema have the same column names.
+
+    # Helper functions for
+    def is_str_type(typ: pa.DataType) -> bool:
+        return cast(bool, typ == pa.string() or typ == pa.large_string())
+
+    def is_str_col(field: pa.Field) -> bool:
+        return is_str_type(field.type)
+
+    def is_str_cat_col(field: pa.Field) -> bool:
+        if not pa.types.is_dictionary(field.type):
+            return False
+        return is_str_type(field.type.value_type)
+
+    # Make a quick check of the old and new schemas to see if any columns need
+    # changing between plain string and categorical-of-string.  We're about to
+    # duplicate the new data -- and we must, since pyarrow.Table is immutable --
+    # so let's only do that if we need to.
+    any_to_change = False
+    for name in new_schema.names:
+        if name not in old_schema.names:
+            continue
+        if is_str_col(old_schema.field(name)) and is_str_cat_col(
+            new_schema.field(name)
+        ):
+            any_to_change = True
+            break
+        if is_str_cat_col(old_schema.field(name)) and is_str_col(
+            new_schema.field(name)
+        ):
+            any_to_change = True
+            break
+
+    if any_to_change:
+        fields_dict = {}
+        for name in new_schema.names:
+            if name not in old_schema.names:
+                continue
+            column = arrow_table.column(name)
+            old_info = old_schema.field(name)
+            new_info = new_schema.field(name)
+            if is_str_col(old_info) and is_str_cat_col(new_info):
+                # Convert from categorical-of-string to plain string.
+                column = column.to_pylist()
+            elif is_str_cat_col(old_info) and is_str_col(new_info):
+                # Convert from plain string to categorical-of-string.  Note:
+                # libtiledbsoma already merges the enum mappings, e.g if the
+                # storage has red, yellow, & green, but our new data has some
+                # yellow, green, and orange.
+                column = pa.array(
+                    column.to_pylist(),
+                    pa.dictionary(
+                        index_type=old_info.type.index_type,
+                        value_type=old_info.type.value_type,
+                        ordered=old_info.type.ordered,
+                    ),
+                )
+            fields_dict[name] = column
+        arrow_table = pa.Table.from_pydict(fields_dict)
+
+    return arrow_table
+
+
 def _extract_new_values_for_append(
     df_uri: str,
     arrow_table: pa.Table,
@@ -1098,15 +1224,10 @@ def _extract_new_values_for_append(
         with _factory.open(
             df_uri, "r", soma_type=DataFrame, context=context
         ) as previous_soma_dataframe:
-            previous_table = previous_soma_dataframe.read().concat()
-            previous_df = previous_table.to_pandas()
-            previous_join_ids = set(
-                int(e) for e in get_dataframe_values(previous_df, SOMA_JOINID)
+            return _extract_new_values_for_append_aux(
+                previous_soma_dataframe, arrow_table
             )
-            mask = [
-                e.as_py() not in previous_join_ids for e in arrow_table[SOMA_JOINID]
-            ]
-            return arrow_table.filter(mask)
+
     except DoesNotExistError:
         return arrow_table
 
diff --git a/apis/python/src/tiledbsoma/io/spatial/outgest.py b/apis/python/src/tiledbsoma/io/spatial/outgest.py
index cc3fb6e423..7b23ccd7a1 100644
--- a/apis/python/src/tiledbsoma/io/spatial/outgest.py
+++ b/apis/python/src/tiledbsoma/io/spatial/outgest.py
@@ -13,6 +13,13 @@
 except ImportError as err:
     warnings.warn("Experimental spatial outgestor requires the spatialdata package.")
     raise err
+
+try:
+    import dask.dataframe as dd
+except ImportError as err:
+    warnings.warn("Experimental spatial outgestor requires the dask package.")
+    raise err
+
 try:
     import geopandas as gpd
 except ImportError as err:
@@ -96,15 +103,18 @@ def _transform_to_spatial_data(
 def to_spatial_data_points(
     points: PointCloudDataFrame,
     *,
+    key: str,
     scene_id: str,
     scene_dim_map: Dict[str, str],
-    transform: somacore.CoordinateTransform,
+    transform: Optional[somacore.CoordinateTransform],
     soma_joinid_name: str,
-) -> pd.DataFrame:
+) -> dd.DataFrame:
     """Export a :class:`PointCloudDataFrame` to a :class:`spatialdata.ShapesModel.
 
     Args:
         points: The point cloud data frame to convert to SpatialData shapes.
+        key: Key for the item in the SpatialData object. Used to set a transformation
+            to the item itself if no scene transform is provided.
         scene_id: The ID of the scene this point cloud dataframe is from.
         scene_dim_map: A mapping from the axis names of the scene to the corresponding
             SpatialData dimension names.
@@ -119,31 +129,36 @@ def to_spatial_data_points(
 
     # Create the SpatialData transform from the points to the Scene (inverse of the
     # transform SOMA stores).
-    transforms = {
-        scene_id: _transform_to_spatial_data(
-            transform.inverse_transform(), points_dim_map, scene_dim_map
-        )
-    }
+    if transform is None:
+        transforms = {key: sd.transformations.Identity()}
+    else:
+        transforms = {
+            scene_id: _transform_to_spatial_data(
+                transform.inverse_transform(), points_dim_map, scene_dim_map
+            )
+        }
 
     # Read the pandas dataframe, rename SOMA_JOINID, add metadata, and return.
     df: pd.DataFrame = points.read().concat().to_pandas()
     df.rename(columns={SOMA_JOINID: soma_joinid_name}, inplace=True)
-    df.attrs["transform"] = transforms
-    return df
+    return sd.models.PointsModel.parse(df, transformations=transforms)
 
 
 def to_spatial_data_shapes(
     points: PointCloudDataFrame,
     *,
+    key: str,
     scene_id: str,
     scene_dim_map: Dict[str, str],
-    transform: somacore.CoordinateTransform,
+    transform: Optional[somacore.CoordinateTransform],
     soma_joinid_name: str,
 ) -> gpd.GeoDataFrame:
     """Export a :class:`PointCloudDataFrame` to a :class:`spatialdata.ShapesModel.
 
     Args:
         points: The point cloud data frame to convert to SpatialData shapes.
+        key: Key for the item in the SpatialData object. Used to set a transformation
+            to the item itself if no scene transform is provided.
         scene_id: The ID of the scene this point cloud dataframe is from.
         scene_dim_map: A mapping from the axis names of the scene to the corresponding
             SpatialData dimension names.
@@ -176,11 +191,14 @@ def to_spatial_data_shapes(
 
     # Create the SpatialData transform from the points to the Scene (inverse of the
     # transform SOMA stores).
-    transforms = {
-        scene_id: _transform_to_spatial_data(
-            transform.inverse_transform(), points_dim_map, scene_dim_map
-        )
-    }
+    if transform is None:
+        transforms = {key: sd.transformations.Identity()}
+    else:
+        transforms = {
+            scene_id: _transform_to_spatial_data(
+                transform.inverse_transform(), points_dim_map, scene_dim_map
+            )
+        }
 
     data = points.read().concat().to_pandas()
     data.rename(columns={SOMA_JOINID: soma_joinid_name}, inplace=True)
@@ -204,12 +222,24 @@ def to_spatial_data_image(
     image: MultiscaleImage,
     level: Optional[Union[str, int]] = None,
     *,
+    key: str,
     scene_id: str,
     scene_dim_map: Dict[str, str],
-    transform: somacore.CoordinateTransform,
+    transform: Optional[somacore.CoordinateTransform],
 ) -> "DataArray":
     """Export a level of a :class:`MultiscaleImage` to a
     :class:`spatialdata.Image2DModel` or :class:`spatialdata.Image3DModel`.
+
+    Args:
+        image: The multiscale image to convert to a level from to a SpatialData image.
+        level: The level of the multiscale image to convert.
+        key: Key for the item in the SpatialData object. Used to set a transformation
+            to the item itself if no scene transform is provided.
+        scene_id: The ID of the scene this multiscale image is from.
+        scene_dim_map: A mapping from the axis names of the scene to the corresponding
+            SpatialData dimension names.
+        transform: The transformation from the coordinate space of the scene this
+            multiscale image is in to the coordinate space of the image itself.
     """
     if not image.has_channel_axis:
         raise NotImplementedError(
@@ -237,29 +267,37 @@ def to_spatial_data_image(
         level = 0
     level_uri = image.level_uri(level)
 
-    # Get the transformtion from the image level to the scene:
-    # If the result is a single scale transform (or identity transform), output a
-    # single transformation. Otherwise, convert to a SpatialData sequence of
-    # transformations.
-    inv_transform = transform.inverse_transform()
-    scale_transform = image.get_transform_from_level(level)
-    if isinstance(transform, somacore.ScaleTransform) or isinstance(
-        scale_transform, somacore.IdentityTransform
-    ):
-        # inv_transform @ scale_transform -> applies scale_transform first
+    if transform is None:
+        # Get the transformation from the image level to the highest resolution of the multiscale image.
+        scale_transform = image.get_transform_from_level(level)
         sd_transform = _transform_to_spatial_data(
-            inv_transform @ scale_transform, image_dim_map, scene_dim_map
-        )
-    else:
-        sd_transform1 = _transform_to_spatial_data(
             scale_transform, image_dim_map, image_dim_map
         )
-        sd_transform2 = _transform_to_spatial_data(
-            inv_transform, image_dim_map, scene_dim_map
-        )
-        # Sequence([sd_transform1, sd_transform2]) -> applies sd_transform1 first
-        sd_transform = sd.transformations.Sequence([sd_transform1, sd_transform2])
-    transformations = {scene_id: sd_transform}
+        transformations = {key: sd_transform}
+    else:
+        # Get the transformation from the image level to the scene:
+        # If the result is a single scale transform (or identity transform), output a
+        # single transformation. Otherwise, convert to a SpatialData sequence of
+        # transformations.
+        inv_transform = transform.inverse_transform()
+        scale_transform = image.get_transform_from_level(level)
+        if isinstance(transform, somacore.ScaleTransform) or isinstance(
+            scale_transform, somacore.IdentityTransform
+        ):
+            # inv_transform @ scale_transform -> applies scale_transform first
+            sd_transform = _transform_to_spatial_data(
+                inv_transform @ scale_transform, image_dim_map, scene_dim_map
+            )
+        else:
+            sd_transform1 = _transform_to_spatial_data(
+                scale_transform, image_dim_map, image_dim_map
+            )
+            sd_transform2 = _transform_to_spatial_data(
+                inv_transform, image_dim_map, scene_dim_map
+            )
+            # Sequence([sd_transform1, sd_transform2]) -> applies sd_transform1 first
+            sd_transform = sd.transformations.Sequence([sd_transform1, sd_transform2])
+        transformations = {scene_id: sd_transform}
 
     # Return array accessor as a dask array.
     return dense_nd_array_to_data_array(
@@ -273,11 +311,23 @@ def to_spatial_data_image(
 def to_spatial_data_multiscale_image(
     image: MultiscaleImage,
     *,
+    key: str,
     scene_id: str,
     scene_dim_map: Dict[str, str],
-    transform: somacore.CoordinateTransform,
+    transform: Optional[somacore.CoordinateTransform],
 ) -> "DataTree":
-    """Export a MultiscaleImage to a DataTree."""
+    """Export a MultiscaleImage to a DataTree.
+
+    Args:
+        image: The multiscale image to convert to a SpatialData image.
+        key: Key for the item in the SpatialData object. Used to set a transformation
+            to the item itself if no scene transform is provided.
+        scene_id: The ID of the scene this multiscale image is from.
+        scene_dim_map: A mapping from the axis names of the scene to the corresponding
+            SpatialData dimension names.
+        transform: The transformation from the coordinate space of the scene this
+            multiscale image is in to the coordinate space of the image itself.
+    """
 
     # Check for channel axis.
     if not image.has_channel_axis:
@@ -296,47 +346,60 @@ def to_spatial_data_multiscale_image(
         orig_axis_names, image.data_axis_order
     )
 
-    # Get the transformtion from the image level to the scene:
-    # If the result is a single scale transform (or identity transform), output a
-    # single transformation. Otherwise, convert to a SpatialData sequence of
-    # transformations.
-    inv_transform = transform.inverse_transform()
-    if isinstance(transform, somacore.ScaleTransform):
-        # inv_transform @ scale_transform -> applies scale_transform first
+    if transform is None:
         spatial_data_transformations = tuple(
             _transform_to_spatial_data(
-                inv_transform @ image.get_transform_from_level(level),
+                image.get_transform_from_level(level),
+                image_dim_map,
                 image_dim_map,
-                scene_dim_map,
             )
             for level in range(image.level_count)
         )
-
     else:
-        sd_scale_transforms = tuple(
-            _transform_to_spatial_data(
-                image.get_transform_from_level(level), image_dim_map, image_dim_map
+        # Get the transformtion from the image level to the scene:
+        # If the result is a single scale transform (or identity transform), output a
+        # single transformation. Otherwise, convert to a SpatialData sequence of
+        # transformations.
+        inv_transform = transform.inverse_transform()
+        if isinstance(transform, somacore.ScaleTransform):
+            # inv_transform @ scale_transform -> applies scale_transform first
+            spatial_data_transformations = tuple(
+                _transform_to_spatial_data(
+                    inv_transform @ image.get_transform_from_level(level),
+                    image_dim_map,
+                    scene_dim_map,
+                )
+                for level in range(image.level_count)
+            )
+        else:
+            sd_scale_transforms = tuple(
+                _transform_to_spatial_data(
+                    image.get_transform_from_level(level), image_dim_map, image_dim_map
+                )
+                for level in range(1, image.level_count)
+            )
+            sd_inv_transform = _transform_to_spatial_data(
+                inv_transform, image_dim_map, scene_dim_map
             )
-            for level in range(1, image.level_count)
-        )
-        sd_inv_transform = _transform_to_spatial_data(
-            inv_transform, image_dim_map, scene_dim_map
-        )
 
-        # First level transform is always the identity, so just directly use
-        # inv_transform. For remaining transformations,
-        # Sequence([sd_transform1, sd_transform2]) -> applies sd_transform1 first
-        spatial_data_transformations = (sd_inv_transform,) + tuple(
-            sd.transformations.Sequence([scale_transform, sd_inv_transform])
-            for scale_transform in sd_scale_transforms
-        )
+            # First level transform is always the identity, so just directly use
+            # inv_transform. For remaining transformations,
+            # Sequence([sd_transform1, sd_transform2]) -> applies sd_transform1 first
+            spatial_data_transformations = (sd_inv_transform,) + tuple(
+                sd.transformations.Sequence([scale_transform, sd_inv_transform])
+                for scale_transform in sd_scale_transforms
+            )
 
     # Create a sequence of resolution level.
     image_data_arrays = tuple(
         dense_nd_array_to_data_array(
             uri=image.level_uri(index),
             dim_names=new_axis_names,
-            attrs={"transform": {scene_id: spatial_data_transformations[index]}},
+            attrs=(
+                {"transform": {key: spatial_data_transformations[index]}}
+                if transform is None
+                else {"transform": {scene_id: spatial_data_transformations[index]}}
+            ),
             context=image.context,
         )
         for index, (soma_name, val) in enumerate(image.levels().items())
diff --git a/apis/python/src/tiledbsoma/managed_query.cc b/apis/python/src/tiledbsoma/managed_query.cc
index 7de08619bc..e3865f5237 100644
--- a/apis/python/src/tiledbsoma/managed_query.cc
+++ b/apis/python/src/tiledbsoma/managed_query.cc
@@ -117,7 +117,10 @@ void load_managed_query(py::module& m) {
             "names"_a,
             "if_not_empty"_a = false)
 
-        .def("submit_read", &ManagedQuery::submit_read)
+        .def(
+            "submit_read",
+            &ManagedQuery::submit_read,
+            py::call_guard<py::gil_scoped_release>())
         .def(
             "results",
             [](ManagedQuery& mq) -> std::optional<py::object> {
@@ -133,6 +136,27 @@ void load_managed_query(py::module& m) {
                 }
             })
 
+        .def(
+            "next",
+            [](ManagedQuery& mq) -> std::optional<py::object> {
+                // Release python GIL before reading data
+                py::gil_scoped_release release;
+                std::optional<std::shared_ptr<ArrayBuffers>> tbl;
+                try {
+                    tbl = mq.read_next();
+                    // Acquire python GIL before accessing python objects
+                } catch (const std::exception& e) {
+                    throw TileDBSOMAError(e.what());
+                }
+                py::gil_scoped_acquire acquire;
+
+                if (!tbl) {
+                    throw py::stop_iteration();
+                }
+
+                return to_table(tbl);
+            })
+
         .def(
             "set_array_data",
             [](ManagedQuery& mq, py::handle py_batch) {
@@ -143,6 +167,7 @@ void load_managed_query(py::module& m) {
                 py_batch.attr("_export_to_c")(
                     arrow_array_ptr, arrow_schema_ptr);
 
+                py::gil_scoped_release release;
                 try {
                     mq.set_array_data(
                         std::make_unique<ArrowSchema>(arrow_schema),
@@ -150,6 +175,7 @@ void load_managed_query(py::module& m) {
                 } catch (const std::exception& e) {
                     TPY_ERROR_LOC(e.what());
                 }
+                py::gil_scoped_acquire acquire;
 
                 arrow_schema.release(&arrow_schema);
                 arrow_array.release(&arrow_array);
@@ -158,17 +184,21 @@ void load_managed_query(py::module& m) {
             "set_soma_data",
             [](ManagedQuery& mq, py::array data) {
                 py::buffer_info data_info = data.request();
+
+                py::gil_scoped_release release;
                 mq.setup_write_column(
                     "soma_data",
                     data.size(),
                     (const void*)data_info.ptr,
                     static_cast<uint64_t*>(nullptr),
                     static_cast<uint8_t*>(nullptr));
+                py::gil_scoped_acquire acquire;
             })
         .def(
             "submit_write",
             &ManagedQuery::submit_write,
-            "sort_coords"_a = false)
+            "sort_coords"_a = false,
+            py::call_guard<py::gil_scoped_release>())
 
         .def("reset", &ManagedQuery::reset)
         .def("close", &ManagedQuery::close)
diff --git a/apis/python/src/tiledbsoma/options/_soma_tiledb_context.py b/apis/python/src/tiledbsoma/options/_soma_tiledb_context.py
index 163578080f..306352316e 100644
--- a/apis/python/src/tiledbsoma/options/_soma_tiledb_context.py
+++ b/apis/python/src/tiledbsoma/options/_soma_tiledb_context.py
@@ -10,7 +10,7 @@
 import threading
 import time
 from concurrent.futures import ThreadPoolExecutor
-from typing import Any, Dict, Literal, Mapping, Optional, Union
+from typing import Any, Dict, Literal, Mapping
 
 from somacore import ContextBase
 from typing_extensions import Self
@@ -19,16 +19,36 @@
 from .._types import OpenTimestamp
 from .._util import ms_to_datetime, to_timestamp_ms
 
+try:
+    from tiledb import Ctx as TileDBCtx
+
+    TILEDB_EXISTS = True
+except ModuleNotFoundError:
+    # If we set this to None, then the type hint TileDBCtx | None in the code
+    # below will error out with:
+    #   TypeError: unsupported operand type(s) for |: 'NoneType' and 'NoneType'
+    # We prevent this by using Any as a catch-all type. This allows the type hint
+    # to remain compatible even when TileDB is not installed, without causing any
+    # issues at runtime or with type checking tools
+    TileDBCtx = Any
+    TILEDB_EXISTS = False
+
+
+def _check_tiledb_ctx() -> None:
+    if not TILEDB_EXISTS:
+        raise ModuleNotFoundError(
+            "The 'tiledb' module is required to access 'tiledb_ctx' but is "
+            "not installed."
+        )
+
 
-def _default_config(
-    override: Mapping[str, Union[str, float]]
-) -> Dict[str, Union[str, float]]:
+def _default_config(override: Mapping[str, str | float]) -> Dict[str, str | float]:
     """Returns a fresh dictionary with TileDB config values.
 
     These should be reasonable defaults that can be used out-of-the-box.
     ``override`` does exactly what it says: overrides default entries.
     """
-    cfg: Dict[str, Union[str, float]] = {
+    cfg: Dict[str, str | float] = {
         "sm.mem.reader.sparse_global_order.ratio_array_data": 0.3
     }
     cfg.update(override)
@@ -41,7 +61,7 @@ def _default_global_native_context() -> clib.SOMAContext:
     return clib.SOMAContext({k: str(v) for k, v in _default_config({}).items()})
 
 
-def _maybe_timestamp_ms(input: Optional[OpenTimestamp]) -> Optional[int]:
+def _maybe_timestamp_ms(input: OpenTimestamp | None) -> int | None:
     if input is None:
         return None
     return to_timestamp_ms(input)
@@ -66,26 +86,35 @@ class SOMATileDBContext(ContextBase):
 
     def __init__(
         self,
-        tiledb_config: Optional[Dict[str, Union[str, float]]] = None,
-        timestamp: Optional[OpenTimestamp] = None,
-        threadpool: Optional[ThreadPoolExecutor] = None,
+        tiledb_config: Dict[str, str | float] | None = None,
+        tiledb_ctx: TileDBCtx | None = None,
+        timestamp: OpenTimestamp | None = None,
+        threadpool: ThreadPoolExecutor | None = None,
     ) -> None:
         """Initializes a new SOMATileDBContext.
 
-        Either ``tiledb_ctx`` or ``tiledb_config`` may be provided, or both may
-        be left at their default. If neither are provided, this will use a
-        single shared :class:`tiledb.Ctx` instantiated upon first use.
-        If ``tiledb_ctx`` is provided, that exact :class:`tiledb.Ctx` is used.
+        Either ``tiledb_config`` or ``tiledb_ctx`` may be provided, or both may
+        be left at their default.
+
         If a ``tiledb_config`` is provided (in the form of a ``dict``),
-        it is used to construct a new ``Ctx``.
+        it is used to construct a new ``SOMAContext``.
 
-        Args:
-            tiledb_ctx: An existing TileDB Context for use as the TileDB Context
-                in this configuration.
+        If `a `tiledb_ctx`` is provided, then it uses the configuration options
+        to construct a new ``SOMAContext``. Note that ``SOMAContext`` will create
+        a new ``Context`` object. The ``tiledb_ctx`` option is only available if
+        the `tiledb` module is installed; otherwise, it will throw a ````ModuleNotFoundError``.
 
+        If neither are provided, this will use a single shared :class:`SOMAContext`
+        instantiated upon first use.
+
+        Args:
             tiledb_config: A set of TileDB configuration options to use,
                 overriding the default configuration.
 
+            tiledb_ctx: A TileDB Context where the set of TileDB
+                configuration options are used to override the default
+                configuration.
+
             timestamp: The default timestamp for operations on SOMA objects,
                 provided either as a ``datetime.datetime`` or a number of
                 milliseconds since the Unix epoch.
@@ -117,22 +146,42 @@ def __init__(
                 provided, a new ThreadPoolExecutor will be created with
                 default settings.
         """
+        if tiledb_ctx is not None and tiledb_config is not None:
+            raise ValueError(
+                "only one of tiledb_config or tiledb_ctx"
+                " may be set when constructing a SOMATileDBContext"
+            )
+
+        # A TileDB Context may only be passed if tiledb is installed
+        if tiledb_ctx is not None:
+            _check_tiledb_ctx()
+
         self._lock = threading.Lock()
         """A lock to ensure single initialization of ``_tiledb_ctx``."""
-        self._initial_config: Optional[Dict[str, Union[str, float]]] = (
+
+        self._initial_config = (
             None if tiledb_config is None else _default_config(tiledb_config)
         )
+        """A dictionary of options to override the default TileDB config.
+
+        This includes both the user-provided options and the default options
+        that we provide to TileDB. If this is unset, then either we were
+        provided with a TileDB Ctx, or we need to use The Default Global Ctx.
+        """
+
+        self._tiledb_ctx = tiledb_ctx
+        """The TileDB context passed in by the user. None if not provided."""
 
-        """The TileDB context to use, either provided or lazily constructed."""
         self._timestamp_ms = _maybe_timestamp_ms(timestamp)
 
         self.threadpool = threadpool or ThreadPoolExecutor()
         """User specified threadpool. If None, we'll instantiate one ourselves."""
-        self._native_context: Optional[clib.SOMAContext] = None
+
+        self._native_context: clib.SOMAContext | None = None
         """Lazily construct clib.SOMAContext."""
 
     @property
-    def timestamp_ms(self) -> Optional[int]:
+    def timestamp_ms(self) -> int | None:
         """
         The default timestamp for SOMA operations, as milliseconds since
         the Unix epoch, or ``None`` if not provided.
@@ -140,7 +189,7 @@ def timestamp_ms(self) -> Optional[int]:
         return self._timestamp_ms
 
     @property
-    def timestamp(self) -> Optional[datetime.datetime]:
+    def timestamp(self) -> datetime.datetime | None:
         """
         The default timestamp for SOMA operations, or ``None`` if not provided.
         """
@@ -153,17 +202,53 @@ def native_context(self) -> clib.SOMAContext:
         """The C++ SOMAContext for this SOMA context."""
         with self._lock:
             if self._native_context is None:
-                if self._initial_config is None:
-                    self._native_context = _default_global_native_context()
-                else:
+                # SOMAContext needs to be lazily created
+                if self._initial_config is not None:
+                    # The user passed in a tiledb_config
                     cfg = self._internal_tiledb_config()
                     self._native_context = clib.SOMAContext(
                         {k: str(v) for k, v in cfg.items()}
                     )
-            return self._native_context
+                elif self._tiledb_ctx is not None:
+                    # The user passed in a tiledb_ctx; if tiledb is not installed
+                    # it should be impossible to enter into this block because
+                    # we already check that in the constructor
+                    assert TileDBCtx is not None
+                    cfg = self._tiledb_ctx.config().dict()
+                    self._native_context = clib.SOMAContext(
+                        {k: str(v) for k, v in cfg.items()}
+                    )
+                else:
+                    # The user did not provide settings so create a default
+                    self._native_context = _default_global_native_context()
+
+        # else SOMAContext already exists
+        return self._native_context
 
     @property
-    def tiledb_config(self) -> Dict[str, Union[str, float]]:
+    def tiledb_ctx(self) -> TileDBCtx | None:
+        """The TileDB-Py Context passed in to create the ``SOMATileDBContext``.
+
+        This accessor is only available if tiledb is installed. If
+        ``SOMATileDBContext`` was constructed with a ``tiledb_ctx``, this returns
+        the ``tiledb.Ctx`` that was passed in. If it was not constructed with
+        a ``tiledb_ctx`` (meaning, either using `tiledb_config` or with default
+        settings), then this will return ``None``.
+
+        Note that when constructing ``SOMATileDBContext`` with ``tiledb_ctx``,
+        internally it uses ``tiledb_ctx.config().dict()`` and constructs a new
+        ``tiledb::Context`` object from the configuration options. Meaning, the
+        ``Context`` passed-in as ``tiledb_ctx`` is NOT the same ``Context`` held
+        by ``SOMATileDBContext``.
+
+        If `tiledb` is not installed, this accessor throws a ``ModuleNotFoundError``
+        error.
+        """
+        _check_tiledb_ctx()
+        return self._tiledb_ctx
+
+    @property
+    def tiledb_config(self) -> Dict[str, str | float]:
         """The TileDB configuration dictionary for this SOMA context.
 
         If this ``SOMATileDBContext`` already has a ``tiledb_ctx``, this will
@@ -177,7 +262,7 @@ def tiledb_config(self) -> Dict[str, Union[str, float]]:
         with self._lock:
             return self._internal_tiledb_config()
 
-    def _internal_tiledb_config(self) -> Dict[str, Union[str, float]]:
+    def _internal_tiledb_config(self) -> Dict[str, str | float]:
         """Internal function for getting the TileDB Config.
 
         Returns a new dict with the contents. Caller must hold ``_lock``.
@@ -186,6 +271,11 @@ def _internal_tiledb_config(self) -> Dict[str, Union[str, float]]:
         if self._native_context is not None:
             return dict(self._native_context.config())
 
+        # The user passed in a TileDB Context. Return its actual config.
+        if self._tiledb_ctx is not None:
+            _check_tiledb_ctx()
+            return dict(self._tiledb_ctx.config())
+
         # Our context has not yet been built.
         # We return what will be passed into the context.
         return (
@@ -197,19 +287,20 @@ def _internal_tiledb_config(self) -> Dict[str, Union[str, float]]:
     def replace(
         self,
         *,
-        tiledb_config: Optional[Dict[str, Any]] = None,
-        timestamp: Union[None, OpenTimestamp, _Unset] = _UNSET,
-        threadpool: Union[None, ThreadPoolExecutor, _Unset] = _UNSET,
+        tiledb_config: Dict[str, Any] | None = None,
+        tiledb_ctx: TileDBCtx | None = None,
+        timestamp: OpenTimestamp | _Unset | None = _UNSET,
+        threadpool: ThreadPoolExecutor | _Unset | None = _UNSET,
     ) -> Self:
         """Create a copy of the context, merging changes.
 
         Args:
             tiledb_config:
-                A dictionary of parameters for `tiledb.Config() <https://tiledb-inc-tiledb.readthedocs-hosted.com/projects/tiledb-py/en/stable/python-api.html#config>`_.
+                A dictionary of parameters representing `tiledb.Config() <https://tiledb-inc-tiledb.readthedocs-hosted.com/projects/tiledb-py/en/stable/python-api.html#config>`_.
                 To remove a parameter from the existing config, provide ``None``
                 as the value.
             tiledb_ctx:
-                A TileDB Context to replace the current context with.
+                A TileDB Context to pull the configuration options from.
             timestamp:
                 A timestamp to replace the current timestamp with.
                 Explicitly passing ``None`` will remove the timestamp.
@@ -230,10 +321,18 @@ def replace(
         """
         with self._lock:
             if tiledb_config is not None:
+                if tiledb_ctx:
+                    raise ValueError(
+                        "Either tiledb_config or tiledb_ctx may be provided"
+                        " to replace(), but not both."
+                    )
                 new_config = self._internal_tiledb_config()
                 new_config.update(tiledb_config)
                 tiledb_config = {k: v for (k, v) in new_config.items() if v is not None}
 
+            if tiledb_ctx is not None:
+                _check_tiledb_ctx()
+
             if timestamp == _UNSET:
                 # Keep the existing timestamp if not overridden.
                 timestamp = self._timestamp_ms
@@ -244,11 +343,12 @@ def replace(
         assert timestamp is None or isinstance(timestamp, (datetime.datetime, int))
         return type(self)(
             tiledb_config=tiledb_config,
+            tiledb_ctx=tiledb_ctx,
             timestamp=timestamp,
             threadpool=threadpool,
         )
 
-    def _open_timestamp_ms(self, in_timestamp: Optional[OpenTimestamp]) -> int:
+    def _open_timestamp_ms(self, in_timestamp: OpenTimestamp | None) -> int:
         """Returns the real timestamp that should be used to open an object."""
         if in_timestamp is not None:
             return to_timestamp_ms(in_timestamp)
@@ -268,6 +368,11 @@ def _validate_soma_tiledb_context(context: Any) -> SOMATileDBContext:
     if context is None:
         return SOMATileDBContext()
 
+    if TILEDB_EXISTS and isinstance(context, TileDBCtx):
+        raise TypeError(
+            "context is a tiledb.Ctx, not a SOMATileDBContext -- please wrap it in tiledbsoma.SOMATileDBContext(...)"
+        )
+
     if not isinstance(context, SOMATileDBContext):
         raise TypeError("context is not a SOMATileDBContext")
 
diff --git a/apis/python/src/tiledbsoma/soma_collection.cc b/apis/python/src/tiledbsoma/soma_collection.cc
index 6c3dd54e8e..5836ca54bb 100644
--- a/apis/python/src/tiledbsoma/soma_collection.cc
+++ b/apis/python/src/tiledbsoma/soma_collection.cc
@@ -76,8 +76,22 @@ void load_soma_collection(py::module& m) {
     py::class_<SOMAMeasurement, SOMACollection, SOMAGroup, SOMAObject>(
         m, "SOMAMeasurement");
 
-    py::class_<SOMAScene, SOMACollection, SOMAGroup, SOMAObject>(
-        m, "SOMAScene");
+    py::class_<SOMAScene, SOMACollection, SOMAGroup, SOMAObject>(m, "SOMAScene")
+        .def_static(
+            "create",
+            [](std::shared_ptr<SOMAContext> ctx,
+               std::string_view uri,
+               std::optional<TimestampRange> timestamp) {
+                try {
+                    SOMAScene::create(uri, ctx, timestamp);
+                } catch (const std::exception& e) {
+                    TPY_ERROR_LOC(e.what());
+                }
+            },
+            py::kw_only(),
+            "ctx"_a,
+            "uri"_a,
+            "timestamp"_a = py::none());
 
     py::class_<SOMAMultiscaleImage, SOMACollection, SOMAGroup, SOMAObject>(
         m, "SOMAMultiscaleImage");
diff --git a/apis/python/src/tiledbsoma/soma_context.cc b/apis/python/src/tiledbsoma/soma_context.cc
index f573d9eeb8..babee68640 100644
--- a/apis/python/src/tiledbsoma/soma_context.cc
+++ b/apis/python/src/tiledbsoma/soma_context.cc
@@ -51,6 +51,7 @@ void load_soma_context(py::module& m) {
     py::class_<SOMAContext, std::shared_ptr<SOMAContext>>(m, "SOMAContext")
         .def(py::init<>())
         .def(py::init<std::map<std::string, std::string>>())
-        .def("config", &SOMAContext::tiledb_config);
+        .def("config", &SOMAContext::tiledb_config)
+        .def("tiledb_ctx", &SOMAContext::tiledb_ctx);
 };
 }  // namespace libtiledbsomacpp
diff --git a/apis/r/NEWS.md b/apis/r/NEWS.md
index 1fb4bbe06d..d6e4b85ffc 100644
--- a/apis/r/NEWS.md
+++ b/apis/r/NEWS.md
@@ -36,6 +36,14 @@
 * Run tests in closer to internal-dependency order by @johnkerl in [#3311](https://github.com/single-cell-data/TileDB-SOMA/pull/3311)
 * Use standard name `zzz.R` for init code by @johnkerl in [#3332](https://github.com/single-cell-data/TileDB-SOMA/pull/3332)
 * Address two more compiler warnings by @johnkerl in [#3320](https://github.com/single-cell-data/TileDB-SOMA/pull/3320)
+* Add support for writing Seurat v5 ragged arrays
+* Update docstrings for `domain` argument to `create` [#3396](https://github.com/single-cell-data/TileDB-SOMA/pull/3396)
+* Vignette for new-shape feature [#3302](https://github.com/single-cell-data/TileDB-SOMA/pull/3302)
+* Fix blockwise iterator + re-indexer to return re-indexed shape instead of full domain
+* Docstring audit for new shape [#3399](https://github.com/single-cell-data/TileDB-SOMA/pull/3399)
+* Add `check_only` support for domain/shape updates [#3400](https://github.com/single-cell-data/TileDB-SOMA/pull/3400)
+* Adjust blockwise + re-indexer to return condensed matrix, not full domain [#3395](https://github.com/single-cell-data/TileDB-SOMA/pull/3395)
+* [c++] Use core 2.27.0.rc5 (#3410) [#3410](https://github.com/single-cell-data/TileDB-SOMA/pull/3410)
 
 # tiledbsoma 1.14.5
 
diff --git a/apis/r/R/BlockwiseIter.R b/apis/r/R/BlockwiseIter.R
index ff90391744..68d6a6ccc1 100644
--- a/apis/r/R/BlockwiseIter.R
+++ b/apis/r/R/BlockwiseIter.R
@@ -346,10 +346,26 @@ BlockwiseSparseReadIter <- R6::R6Class(
     .zero_based = FALSE,
     soma_reader_transform = function(x) {
       tbl <- private$reindex_arrow_table(soma_array_to_arrow_table(x))
+      shape <- private$.shape
+      axis <- as.integer(self$axis)
+      axname <- sprintf("soma_dim_%i", axis)
+      stride <- self$coords[[axname]]$stride
+      # For re-indexed blockwise iterators, shape should reflect the re-indexed
+      # axes; this can generally be the stride of the iterator, but for the end
+      # of each iterator this needs to be the remainder (see ?`%%`).
+      # Note: this only applies to the major axis, minor axes always return
+      # the full domain.
+      if (self$reindexable && shape[axis + 1L] > stride) {
+        shape[axis + 1L] <- if (self$coords[[axname]]$has_next()) {
+          as.numeric(stride)
+        } else {
+          as.numeric(self$coords[[axname]]$end %% stride) + 1L
+        }
+      }
       mat <- arrow_table_to_sparse(
         tbl,
         repr = self$repr,
-        shape = private$.shape,
+        shape = shape,
         zero_based = private$.zero_based
       )
       attr(mat, "coords") <- attr(tbl, "coords", exact = TRUE)
diff --git a/apis/r/R/RcppExports.R b/apis/r/R/RcppExports.R
index e49188b65e..cb06831a9a 100644
--- a/apis/r/R/RcppExports.R
+++ b/apis/r/R/RcppExports.R
@@ -238,20 +238,20 @@ c_schema <- function(uri, ctxxp) {
     .Call(`_tiledbsoma_c_schema`, uri, ctxxp)
 }
 
-resize <- function(uri, new_shape, function_name_for_messages, ctxxp) {
-    invisible(.Call(`_tiledbsoma_resize`, uri, new_shape, function_name_for_messages, ctxxp))
+resize <- function(uri, new_shape, function_name_for_messages, check_only, ctxxp) {
+    .Call(`_tiledbsoma_resize`, uri, new_shape, function_name_for_messages, check_only, ctxxp)
 }
 
 resize_soma_joinid_shape <- function(uri, new_shape, function_name_for_messages, ctxxp) {
     invisible(.Call(`_tiledbsoma_resize_soma_joinid_shape`, uri, new_shape, function_name_for_messages, ctxxp))
 }
 
-tiledbsoma_upgrade_shape <- function(uri, new_shape, function_name_for_messages, ctxxp) {
-    invisible(.Call(`_tiledbsoma_tiledbsoma_upgrade_shape`, uri, new_shape, function_name_for_messages, ctxxp))
+tiledbsoma_upgrade_shape <- function(uri, new_shape, function_name_for_messages, check_only, ctxxp) {
+    .Call(`_tiledbsoma_tiledbsoma_upgrade_shape`, uri, new_shape, function_name_for_messages, check_only, ctxxp)
 }
 
-upgrade_or_change_domain <- function(uri, is_change_domain, nadimap, nadimsp, function_name_for_messages, ctxxp) {
-    invisible(.Call(`_tiledbsoma_upgrade_or_change_domain`, uri, is_change_domain, nadimap, nadimsp, function_name_for_messages, ctxxp))
+upgrade_or_change_domain <- function(uri, is_change_domain, nadimap, nadimsp, function_name_for_messages, check_only, ctxxp) {
+    .Call(`_tiledbsoma_upgrade_or_change_domain`, uri, is_change_domain, nadimap, nadimsp, function_name_for_messages, check_only, ctxxp)
 }
 
 c_update_dataframe_schema <- function(uri, ctxxp, column_names_to_drop, add_cols_types, add_cols_enum_value_types, add_cols_enum_ordered) {
diff --git a/apis/r/R/SOMAArrayBase.R b/apis/r/R/SOMAArrayBase.R
index 9a49db817f..5514a7d915 100644
--- a/apis/r/R/SOMAArrayBase.R
+++ b/apis/r/R/SOMAArrayBase.R
@@ -1,7 +1,7 @@
 #' SOMA Array Base Class
 #'
 #' Adds SOMA-specific functionality to the [`TileDBArray`] class.  (lifecycle:
-#' experimental)
+#' maturing)
 #' @keywords internal
 
 SOMAArrayBase <- R6::R6Class(
diff --git a/apis/r/R/SOMADataFrame.R b/apis/r/R/SOMADataFrame.R
index 761e0d98ff..06b1d23b55 100644
--- a/apis/r/R/SOMADataFrame.R
+++ b/apis/r/R/SOMADataFrame.R
@@ -19,14 +19,22 @@ SOMADataFrame <- R6::R6Class(
     #' @param index_column_names A vector of column names to use as user-defined
     #' index columns.  All named columns must exist in the schema, and at least
     #' one index column name is required.
-    #' @param domain An optional list of 2-element vectors specifying the domain of each index
-    #' column. Each vector should be a pair consisting of the minimum and maximum values storable in
-    #' the index column. For example, if there is a single int64-valued index column, then `domain`
-    #' might be `c(100, 200)` to indicate that values between 100 and 200, inclusive, can be stored
-    #' in that column.  If provided, this list must have the same length as `index_column_names`,
-    #' and the index-column domain will be as specified.  If omitted entirely, or if `NULL` in a given
-    #' dimension, the corresponding index-column domain will use the minimum and maximum possible
-    #' values for the column's datatype.  This makes a `DataFrame` growable.
+    #' @param domain An optional list specifying the domain of each
+    #' index column. Each slot in the list must have its name being the name
+    #' of an index column, and its value being be a length-two vector
+    #' consisting of the minimum and maximum values storable in the index
+    #' column. For example, if there is a single int64-valued index column
+    #' `soma_joinid`, then `domain` might be `list(soma_joinid=c(100, 200))`
+    #' to indicate that values between 100 and 200, inclusive, can be stored
+    #' in that column.  If provided, this sequence must have the same length
+    #' as `index_column_names`, and the index-column domain will be as
+    #' specified.  If omitted entirely, or if `NULL` in a given dimension, the
+    #' corresponding index-column domain will use an empty range, and data writes
+    #' after that will fail with "A range was set outside of the current domain".
+    #' Unless you have a particular reason not to, you should always provide the
+    #' desired `domain` at create time: this is an optional but strongly
+    #' recommended parameter.  See also `change_domain` which allows you to
+    #' expand the domain after create.
     #' @template param-platform-config
     #' @param internal_use_only Character value to signal this is a 'permitted' call,
     #' as `create()` is considered internal and should not be called directly.
@@ -464,9 +472,11 @@ SOMADataFrame <- R6::R6Class(
     #' @param new_domain A named list, keyed by index-column name, with values
     #' being two-element vectors containing the desired lower and upper bounds
     #' for the domain.
-    #' @return No return value
-    tiledbsoma_upgrade_domain = function(new_domain) {
-      # stopifnot("'new_domain' must be CODE ME UP PLZ" = ...
+    #' @param check_only If true, does not apply the operation, but only reports
+    #' whether it would have succeeded.
+    #' @return No return value if `check_only` is `FALSE`. If `check_only` is `TRUE`,
+    #' returns the empty string if no error is detected, else a description of the error.
+    tiledbsoma_upgrade_domain = function(new_domain, check_only = FALSE) {
       # Checking slotwise new shape >= old shape, and <= max_shape, is already
       # done in libtiledbsoma
 
@@ -474,16 +484,23 @@ SOMADataFrame <- R6::R6Class(
         new_domain, "tiledbsoma_upgrade_domain"
       )
 
-      invisible(
-        upgrade_or_change_domain(
-          self$uri,
-          FALSE,
-          pyarrow_domain_table$array,
-          pyarrow_domain_table$schema,
-          .name_of_function(),
-          private$.soma_context
-        )
+      reason_string <- upgrade_or_change_domain(
+        self$uri,
+        FALSE,
+        pyarrow_domain_table$array,
+        pyarrow_domain_table$schema,
+        .name_of_function(),
+        check_only,
+        private$.soma_context
       )
+
+      if (isTRUE(check_only)) {
+        return(reason_string)
+      }
+
+      # The return value from upgrade_or_change_domain without check_only is
+      # always "", or it raises an error trying.
+      return(invisible(NULL))
     },
 
     #' @description Allows you to set the domain of a `SOMADataFrame`, when the
@@ -494,13 +511,15 @@ SOMADataFrame <- R6::R6Class(
     #' whether the operation would succeed if attempted, and a reason why it
     #' would not. The return value from `domain` must be contained within
     #' the requested `new_domain`, and the requested `new_domain` must be
-    #' contained within the return value from `maxdomain`.
+    #' contained within the return value from `maxdomain`. (lifecycle: maturing)
     #' @param new_domain A named list, keyed by index-column name, with values
     #' being two-element vectors containing the desired lower and upper bounds
     #' for the domain.
-    #' @return No return value
-    change_domain = function(new_domain) {
-      # stopifnot("'new_domain' must be CODE ME UP PLZ" = ...
+    #' @param check_only If true, does not apply the operation, but only reports
+    #' whether it would have succeeded.
+    #' @return No return value if `check_only` is `FALSE`. If `check_only` is `TRUE`,
+    #' returns the empty string if no error is detected, else a description of the error.
+    change_domain = function(new_domain, check_only = FALSE) {
       # Checking slotwise new shape >= old shape, and <= max_shape, is already
       # done in libtiledbsoma
 
@@ -508,16 +527,23 @@ SOMADataFrame <- R6::R6Class(
         new_domain, tiledbsoma_upgrade_domain
       )
 
-      invisible(
-        upgrade_or_change_domain(
+      reason_string <- upgrade_or_change_domain(
           self$uri,
           TRUE,
           pyarrow_domain_table$array,
           pyarrow_domain_table$schema,
           .name_of_function(),
+          check_only,
           private$.soma_context
-        )
       )
+
+      if (isTRUE(check_only)) {
+        return(reason_string)
+      }
+
+      # The return value from upgrade_or_change_domain without check_only is
+      # always "", or it raises an error trying.
+      return(invisible(NULL))
     }
   ),
   private = list(
diff --git a/apis/r/R/SOMANDArrayBase.R b/apis/r/R/SOMANDArrayBase.R
index 2f3d5d27a7..f4049a165a 100644
--- a/apis/r/R/SOMANDArrayBase.R
+++ b/apis/r/R/SOMANDArrayBase.R
@@ -14,7 +14,7 @@ SOMANDArrayBase <- R6::R6Class(
   public = list(
 
     #' @description Create a SOMA NDArray named with the URI. (lifecycle:
-    #' experimental)
+    #' maturing)
     #' @param type an [Arrow type][arrow::data-type] defining the type of each
     #' element in the array.
     #' @param shape a vector of integers defining the shape of the array.
@@ -79,7 +79,7 @@ SOMANDArrayBase <- R6::R6Class(
     },
 
     ## needed eg after open() to set (Arrow) type
-    #' @description Sets a cache value for the datatype (lifecycle: experimental)
+    #' @description Sets a cache value for the datatype (lifecycle: maturing)
     #' @param type A character value describing the TileDB data type
     set_data_type = function(type) {
       spdl::debug("[SOMANDArrayBase::set_data_type] caching type {}", type$ToString())
@@ -100,17 +100,29 @@ SOMANDArrayBase <- R6::R6Class(
     #' the current shape in any dimension.  Raises an error if the requested new
     #' shape exceeds `maxshape` in any dimension. Raises an error if the array
     #' doesn't already have a shape: in that case please call
-    #' `tiledbsoma_upgrade_shape`.
+    #' `tiledbsoma_upgrade_shape`. (lifecycle: maturing)
     #' @param new_shape A vector of integerish, of the same length as the array's `ndim`.
-    #' @return No return value
-    resize = function(new_shape) {
+    #' @param check_only If true, does not apply the operation, but only reports
+    #' whether it would have succeeded.
+    #' @return No return value if `check_only` is `FALSE`. If `check_only` is `TRUE`,
+    #' returns the empty string if no error is detected, else a description of the error.
+    resize = function(new_shape, check_only = FALSE) {
       stopifnot(
         "'new_shape' must be a vector of integerish values, of the same length as maxshape" =
           rlang::is_integerish(new_shape, n = self$ndim()) ||
             (bit64::is.integer64(new_shape) && length(new_shape) == self$ndim())
       )
       # Checking slotwise new shape >= old shape, and <= max_shape, is already done in libtiledbsoma
-      resize(self$uri, new_shape, .name_of_function(), private$.soma_context)
+
+      reason_string <- resize(self$uri, new_shape, .name_of_function(), check_only, private$.soma_context)
+
+      if (isTRUE(check_only)) {
+        return(reason_string)
+      }
+      
+      # The return value from resize without check_only is always "", or it
+      # raises an error trying.
+      return(invisible(NULL))
     },
 
     #' @description Allows the array to have a resizeable shape as described in the
@@ -119,16 +131,28 @@ SOMANDArrayBase <- R6::R6Class(
     #' `tiledbsoma_upgrade_shape` and `resize` are very similar: the former must be
     #' call on a pre-1.15 array the first time a shape is set on it; the latter must
     #' be used for subsequent resizes on any array which already has upgraded shape.
+    #' (lifecycle: maturing)
     #' @param shape A vector of integerish, of the same length as the array's `ndim`.
-    #' @return No return value
-    tiledbsoma_upgrade_shape = function(shape) {
+    #' @param check_only If true, does not apply the operation, but only reports
+    #' whether it would have succeeded.
+    #' @return No return value if `check_only` is `FALSE`. If `check_only` is `TRUE`,
+    #' returns the empty string if no error is detected, else a description of the error.
+    tiledbsoma_upgrade_shape = function(shape, check_only = FALSE) {
       stopifnot(
         "'shape' must be a vector of integerish values, of the same length as maxshape" =
           rlang::is_integerish(shape, n = self$ndim()) ||
             (bit64::is.integer64(shape) && length(shape) == self$ndim())
       )
       # Checking slotwise new shape >= old shape, and <= max_shape, is already done in libtiledbsoma
-      tiledbsoma_upgrade_shape(self$uri, shape, .name_of_function(), private$.soma_context)
+
+      reason_string <- tiledbsoma_upgrade_shape(self$uri, shape, .name_of_function(), check_only, private$.soma_context)
+      if (isTRUE(check_only)) {
+        return(reason_string)
+      }
+
+      # The return value from tiledbsoma_upgrade_shape without check_only is
+      # always "", or it raises an error trying.
+      return(invisible(NULL))
     }
   ),
   private = list(
diff --git a/apis/r/R/datasets.R b/apis/r/R/datasets.R
index b09b13a88b..e78a0889d9 100644
--- a/apis/r/R/datasets.R
+++ b/apis/r/R/datasets.R
@@ -50,6 +50,7 @@ extract_dataset <- function(name, dir = tempdir()) {
 
   # Extract tar.gz file to dir
   tarfile <- dir(data_dir, pattern = name, full.names = TRUE)
+  tarfile <- file.path(data_dir, sprintf("%s.tar.gz", name))
   stopifnot("The specified dataset does not exist" = file.exists(tarfile))
 
   dataset_uri <- file.path(dir, name)
diff --git a/apis/r/_pkgdown.yml b/apis/r/_pkgdown.yml
index 9372272e0c..b0ba1f6786 100644
--- a/apis/r/_pkgdown.yml
+++ b/apis/r/_pkgdown.yml
@@ -17,6 +17,8 @@ navbar:
         href: articles/soma-objects.html
       - text: Reading From SOMA Objects
         href: articles/soma-reading.html
+      - text: SOMA shapes
+        href: articles/soma-shapes.html
       - text: Querying a SOMA Experiment
         href: articles/soma-experiment-queries.html
       - text: -------
diff --git a/apis/r/man/SOMAArrayBase.Rd b/apis/r/man/SOMAArrayBase.Rd
index 246111ddac..6edae04ccf 100644
--- a/apis/r/man/SOMAArrayBase.Rd
+++ b/apis/r/man/SOMAArrayBase.Rd
@@ -5,7 +5,7 @@
 \title{SOMA Array Base Class}
 \description{
 Adds SOMA-specific functionality to the \code{\link{TileDBArray}} class.  (lifecycle:
-experimental)
+maturing)
 }
 \keyword{internal}
 \section{Super classes}{
diff --git a/apis/r/man/SOMADataFrame.Rd b/apis/r/man/SOMADataFrame.Rd
index a4ffade826..d85faf0ae4 100644
--- a/apis/r/man/SOMADataFrame.Rd
+++ b/apis/r/man/SOMADataFrame.Rd
@@ -24,7 +24,9 @@ row and is intended to act as a join key for other objects, such as
 \item \href{#method-SOMADataFrame-domain}{\code{SOMADataFrame$domain()}}
 \item \href{#method-SOMADataFrame-maxdomain}{\code{SOMADataFrame$maxdomain()}}
 \item \href{#method-SOMADataFrame-tiledbsoma_has_upgraded_domain}{\code{SOMADataFrame$tiledbsoma_has_upgraded_domain()}}
-\item \href{#method-SOMADataFrame-resize_soma_joinid_shape}{\code{SOMADataFrame$resize_soma_joinid_shape()}}
+\item \href{#method-SOMADataFrame-tiledbsoma_resize_soma_joinid_shape}{\code{SOMADataFrame$tiledbsoma_resize_soma_joinid_shape()}}
+\item \href{#method-SOMADataFrame-tiledbsoma_upgrade_domain}{\code{SOMADataFrame$tiledbsoma_upgrade_domain()}}
+\item \href{#method-SOMADataFrame-change_domain}{\code{SOMADataFrame$change_domain()}}
 \item \href{#method-SOMADataFrame-clone}{\code{SOMADataFrame$clone()}}
 }
 }
@@ -81,14 +83,22 @@ Create (lifecycle: maturing)
 index columns.  All named columns must exist in the schema, and at least
 one index column name is required.}
 
-\item{\code{domain}}{An optional list of 2-element vectors specifying the domain of each index
-column. Each vector should be a pair consisting of the minimum and maximum values storable in
-the index column. For example, if there is a single int64-valued index column, then \code{domain}
-might be \code{c(100, 200)} to indicate that values between 100 and 200, inclusive, can be stored
-in that column.  If provided, this list must have the same length as \code{index_column_names},
-and the index-column domain will be as specified.  If omitted entirely, or if \code{NULL} in a given
-dimension, the corresponding index-column domain will use the minimum and maximum possible
-values for the column's datatype.  This makes a \code{DataFrame} growable.}
+\item{\code{domain}}{An optional list specifying the domain of each
+index column. Each slot in the list must have its name being the name
+of an index column, and its value being be a length-two vector
+consisting of the minimum and maximum values storable in the index
+column. For example, if there is a single int64-valued index column
+\code{soma_joinid}, then \code{domain} might be \code{list(soma_joinid=c(100, 200))}
+to indicate that values between 100 and 200, inclusive, can be stored
+in that column.  If provided, this sequence must have the same length
+as \code{index_column_names}, and the index-column domain will be as
+specified.  If omitted entirely, or if \code{NULL} in a given dimension, the
+corresponding index-column domain will use an empty range, and data writes
+after that will fail with "A range was set outside of the current domain".
+Unless you have a particular reason not to, you should always provide the
+desired \code{domain} at create time: this is an optional but strongly
+recommended parameter.  See also \code{change_domain} which allows you to
+expand the domain after create.}
 
 \item{\code{platform_config}}{A \link[tiledbsoma:PlatformConfig]{platform configuration}
 object}
@@ -139,8 +149,8 @@ column, and optionally filtered.
 \subsection{Arguments}{
 \if{html}{\out{<div class="arguments">}}
 \describe{
-\item{\code{coords}}{Optional named list of indices specifying the rows to read; each (named)
-list element corresponds to a dimension of the same name.}
+\item{\code{coords}}{Optional named list of indices specifying the rows to read;
+each (named) list element corresponds to a dimension of the same name.}
 
 \item{\code{column_names}}{Optional character vector of column names to return.}
 
@@ -280,9 +290,9 @@ Logical
 }
 }
 \if{html}{\out{<hr>}}
-\if{html}{\out{<a id="method-SOMADataFrame-resize_soma_joinid_shape"></a>}}
-\if{latex}{\out{\hypertarget{method-SOMADataFrame-resize_soma_joinid_shape}{}}}
-\subsection{Method \code{resize_soma_joinid_shape()}}{
+\if{html}{\out{<a id="method-SOMADataFrame-tiledbsoma_resize_soma_joinid_shape"></a>}}
+\if{latex}{\out{\hypertarget{method-SOMADataFrame-tiledbsoma_resize_soma_joinid_shape}{}}}
+\subsection{Method \code{tiledbsoma_resize_soma_joinid_shape()}}{
 Increases the shape of the dataframe on the \code{soma_joinid}
 index column, if it indeed is an index column, leaving all other index
 columns as-is. If the \code{soma_joinid} is not an index column, no change is
@@ -292,7 +302,7 @@ domain expansion.  Raises an error if the dataframe doesn't already have a
 domain: in that case please call \code{tiledbsoma_upgrade_domain} (WIP for
 1.15).
 \subsection{Usage}{
-\if{html}{\out{<div class="r">}}\preformatted{SOMADataFrame$resize_soma_joinid_shape(new_shape)}\if{html}{\out{</div>}}
+\if{html}{\out{<div class="r">}}\preformatted{SOMADataFrame$tiledbsoma_resize_soma_joinid_shape(new_shape)}\if{html}{\out{</div>}}
 }
 
 \subsection{Arguments}{
@@ -308,6 +318,65 @@ No return value
 }
 }
 \if{html}{\out{<hr>}}
+\if{html}{\out{<a id="method-SOMADataFrame-tiledbsoma_upgrade_domain"></a>}}
+\if{latex}{\out{\hypertarget{method-SOMADataFrame-tiledbsoma_upgrade_domain}{}}}
+\subsection{Method \code{tiledbsoma_upgrade_domain()}}{
+Allows you to set the domain of a \code{SOMADataFrame}, when the
+\code{SOMADataFrame} does not have a domain set yet.  The argument must be a
+tuple of pairs of low/high values for the desired domain, one pair per
+index column. For string index columns, you must offer the low/high pair
+as \verb{("", "")}, or as \code{NULL}.  If \code{check_only} is \code{True}, returns
+whether the operation would succeed if attempted, and a reason why it
+would not. The domain being requested must be contained within what
+\code{maxdomain} returns.
+\subsection{Usage}{
+\if{html}{\out{<div class="r">}}\preformatted{SOMADataFrame$tiledbsoma_upgrade_domain(new_domain)}\if{html}{\out{</div>}}
+}
+
+\subsection{Arguments}{
+\if{html}{\out{<div class="arguments">}}
+\describe{
+\item{\code{new_domain}}{A named list, keyed by index-column name, with values
+being two-element vectors containing the desired lower and upper bounds
+for the domain.}
+}
+\if{html}{\out{</div>}}
+}
+\subsection{Returns}{
+No return value
+}
+}
+\if{html}{\out{<hr>}}
+\if{html}{\out{<a id="method-SOMADataFrame-change_domain"></a>}}
+\if{latex}{\out{\hypertarget{method-SOMADataFrame-change_domain}{}}}
+\subsection{Method \code{change_domain()}}{
+Allows you to set the domain of a \code{SOMADataFrame}, when the
+\code{SOMADataFrame} already has a domain set yet.  The argument must be a
+tuple of pairs of low/high values for the desired domain, one pair per
+index column. For string index columns, you must offer the low/high pair
+as \verb{("", "")}, or as \code{NULL}.  If \code{check_only} is \code{True}, returns
+whether the operation would succeed if attempted, and a reason why it
+would not. The return value from \code{domain} must be contained within
+the requested \code{new_domain}, and the requested \code{new_domain} must be
+contained within the return value from \code{maxdomain}. (lifecycle: maturing)
+\subsection{Usage}{
+\if{html}{\out{<div class="r">}}\preformatted{SOMADataFrame$change_domain(new_domain)}\if{html}{\out{</div>}}
+}
+
+\subsection{Arguments}{
+\if{html}{\out{<div class="arguments">}}
+\describe{
+\item{\code{new_domain}}{A named list, keyed by index-column name, with values
+being two-element vectors containing the desired lower and upper bounds
+for the domain.}
+}
+\if{html}{\out{</div>}}
+}
+\subsection{Returns}{
+No return value
+}
+}
+\if{html}{\out{<hr>}}
 \if{html}{\out{<a id="method-SOMADataFrame-clone"></a>}}
 \if{latex}{\out{\hypertarget{method-SOMADataFrame-clone}{}}}
 \subsection{Method \code{clone()}}{
diff --git a/apis/r/man/SOMADenseNDArray.Rd b/apis/r/man/SOMADenseNDArray.Rd
index 807c0d5816..5f1d5dce79 100644
--- a/apis/r/man/SOMADenseNDArray.Rd
+++ b/apis/r/man/SOMADenseNDArray.Rd
@@ -9,10 +9,13 @@ offset (zero-based) \code{int64} integer indexing on each dimension with domain
 \verb{[0, maxInt64)}. The \code{SOMADenseNDArray} has a user-defined schema, which
 includes:
 \itemize{
-\item \strong{type}: a \code{primitive} type, expressed as an Arrow type (e.g., \code{int64},
-\code{float32}, etc), indicating the type of data contained within the array
-\item \strong{shape}: the shape of the array, i.e., number and length of each
-dimension
+\item \strong{type}: A \code{primitive} type, expressed as an Arrow type (e.g., \code{int64},
+\code{float32}, etc), indicating the type of data contained within the array.
+\item \strong{shape}: The shape of the array, i.e., number and length of each
+dimension. This is a soft limit which can be increased using \code{resize},
+up to the \code{maxshape}.
+\item \strong{maxshape}: The hard limit up to which \code{shape} may be increased using
+\code{resize}.
 }
 
 All dimensions must have a positive, non-zero length, and there must be 1 or
@@ -81,7 +84,7 @@ Read as an 'arrow::Table' (lifecycle: maturing)
 \if{html}{\out{<div class="r">}}\preformatted{SOMADenseNDArray$read_arrow_table(
   coords = NULL,
   result_order = "auto",
-  log_level = "warn"
+  log_level = "auto"
 )}\if{html}{\out{</div>}}
 }
 
diff --git a/apis/r/man/SOMANDArrayBase.Rd b/apis/r/man/SOMANDArrayBase.Rd
index 66f61baec5..19c0d0c44f 100644
--- a/apis/r/man/SOMANDArrayBase.Rd
+++ b/apis/r/man/SOMANDArrayBase.Rd
@@ -58,7 +58,7 @@ Adds NDArray-specific functionality to the \code{\link{SOMAArrayBase}} class.
 \if{latex}{\out{\hypertarget{method-SOMANDArrayBase-create}{}}}
 \subsection{Method \code{create()}}{
 Create a SOMA NDArray named with the URI. (lifecycle:
-experimental)
+maturing)
 \subsection{Usage}{
 \if{html}{\out{<div class="r">}}\preformatted{SOMANDArrayBase$create(
   type,
@@ -90,7 +90,7 @@ directly.}
 \if{html}{\out{<a id="method-SOMANDArrayBase-set_data_type"></a>}}
 \if{latex}{\out{\hypertarget{method-SOMANDArrayBase-set_data_type}{}}}
 \subsection{Method \code{set_data_type()}}{
-Sets a cache value for the datatype (lifecycle: experimental)
+Sets a cache value for the datatype (lifecycle: maturing)
 \subsection{Usage}{
 \if{html}{\out{<div class="r">}}\preformatted{SOMANDArrayBase$set_data_type(type)}\if{html}{\out{</div>}}
 }
@@ -123,11 +123,12 @@ Logical
 \if{html}{\out{<a id="method-SOMANDArrayBase-resize"></a>}}
 \if{latex}{\out{\hypertarget{method-SOMANDArrayBase-resize}{}}}
 \subsection{Method \code{resize()}}{
-Increases the shape of the array as specfied. Raises an error
-if the new shape is less than the current shape in any dimension. Raises
-an error if the new shape exceeds maxshape in any dimension. Raises an
-error if the array doesn't already have a shape: in that case please call
-tiledbsoma_upgrade_shape.
+Increases the shape of the array as specified, up to the hard
+limit which is \code{maxshape}. Raises an error if the new shape is less than
+the current shape in any dimension.  Raises an error if the requested new
+shape exceeds \code{maxshape} in any dimension. Raises an error if the array
+doesn't already have a shape: in that case please call
+\code{tiledbsoma_upgrade_shape}. (lifecycle: maturing)
 \subsection{Usage}{
 \if{html}{\out{<div class="r">}}\preformatted{SOMANDArrayBase$resize(new_shape)}\if{html}{\out{</div>}}
 }
@@ -149,7 +150,11 @@ No return value
 \subsection{Method \code{tiledbsoma_upgrade_shape()}}{
 Allows the array to have a resizeable shape as described in the
 TileDB-SOMA 1.15 release notes.  Raises an error if the shape exceeds maxshape in any
-dimension. Raises an error if the array already has a shape.
+dimension. Raises an error if the array already has a shape. The methods
+\code{tiledbsoma_upgrade_shape} and \code{resize} are very similar: the former must be
+call on a pre-1.15 array the first time a shape is set on it; the latter must
+be used for subsequent resizes on any array which already has upgraded shape.
+(lifecycle: maturing)
 \subsection{Usage}{
 \if{html}{\out{<div class="r">}}\preformatted{SOMANDArrayBase$tiledbsoma_upgrade_shape(shape)}\if{html}{\out{</div>}}
 }
diff --git a/apis/r/man/SOMASparseNDArray.Rd b/apis/r/man/SOMASparseNDArray.Rd
index 5915ef7fb8..98cb567f9f 100644
--- a/apis/r/man/SOMASparseNDArray.Rd
+++ b/apis/r/man/SOMASparseNDArray.Rd
@@ -8,8 +8,13 @@
 (zero-based) integer indexing on each dimension. The \code{SOMASparseNDArray} has
 a user-defined schema, which includes:
 \itemize{
-\item type - a \code{primitive} type, expressed as an Arrow type (e.g., \code{int64}, \code{float32}, etc)
-\item shape - the shape of the array, i.e., number and length of each dimension
+\item \strong{type}: A \code{primitive} type, expressed as an Arrow type (e.g., \code{int64},
+\code{float32}, etc), indicating the type of data contained within the array.
+\item \strong{shape}: The shape of the array, i.e., number and length of each
+dimension. This is a soft limit which can be increased using \code{resize},
+up to the \code{maxshape}.
+\item \strong{maxshape}: The hard limit up to which \code{shape} may be increased using
+\code{resize}.
 }
 
 All dimensions must have a positive, non-zero length.
diff --git a/apis/r/man/TileDBArray.Rd b/apis/r/man/TileDBArray.Rd
index 533a9557f7..9052706069 100644
--- a/apis/r/man/TileDBArray.Rd
+++ b/apis/r/man/TileDBArray.Rd
@@ -206,8 +206,10 @@ A named list of \code{\link[tiledb:tiledb_dim]{tiledb::tiledb_dim}} objects
 \if{latex}{\out{\hypertarget{method-TileDBArray-shape}{}}}
 \subsection{Method \code{shape()}}{
 Retrieve the shape, i.e. the capacity of each dimension.
-This will not necessarily match the bounds of occupied cells within the
-array.  Rather, it is the bounds outside of which no data may be written.
+Attempted reads and writes outside the \code{shape} will result in a runtime
+error: this is the purpose of \code{shape}.  This will not necessarily match
+the bounds of occupied cells within the array. Using \code{resize}, this may be
+increased up to the hard limit which \code{maxshape} reports.
 (lifecycle: maturing)
 \subsection{Usage}{
 \if{html}{\out{<div class="r">}}\preformatted{TileDBArray$shape()}\if{html}{\out{</div>}}
@@ -221,7 +223,8 @@ A named vector of dimension length (and the same type as the dimension)
 \if{html}{\out{<a id="method-TileDBArray-maxshape"></a>}}
 \if{latex}{\out{\hypertarget{method-TileDBArray-maxshape}{}}}
 \subsection{Method \code{maxshape()}}{
-Retrieve the maxshape, i.e. maximum possible that the
+Retrieve the hard limit up to which the array may be resized
+using the \code{resize} method.
 (lifecycle: maturing)
 \subsection{Usage}{
 \if{html}{\out{<div class="r">}}\preformatted{TileDBArray$maxshape()}\if{html}{\out{</div>}}
diff --git a/apis/r/man/tiledbsoma_empty_query_condition.Rd b/apis/r/man/tiledbsoma_empty_query_condition.Rd
new file mode 100644
index 0000000000..b12551c1d3
--- /dev/null
+++ b/apis/r/man/tiledbsoma_empty_query_condition.Rd
@@ -0,0 +1,18 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/QueryCondition.R
+\name{tiledbsoma_empty_query_condition}
+\alias{tiledbsoma_empty_query_condition}
+\title{Creates a 'tiledbsoma_query_condition' object}
+\usage{
+tiledbsoma_empty_query_condition(somactx)
+}
+\arguments{
+\item{somactx}{(optional) A TileDB Ctx object; if not supplied the default
+context object is retrieved}
+}
+\value{
+A 'tiledbsoma_query_condition' object
+}
+\description{
+Creates a 'tiledbsoma_query_condition' object
+}
diff --git a/apis/r/man/tiledbsoma_query_condition-class.Rd b/apis/r/man/tiledbsoma_query_condition-class.Rd
new file mode 100644
index 0000000000..269456e725
--- /dev/null
+++ b/apis/r/man/tiledbsoma_query_condition-class.Rd
@@ -0,0 +1,18 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/QueryCondition.R
+\docType{class}
+\name{tiledbsoma_query_condition-class}
+\alias{tiledbsoma_query_condition-class}
+\title{An S4 class for a TileDB QueryCondition object}
+\description{
+An S4 class for a TileDB QueryCondition object
+}
+\section{Slots}{
+
+\describe{
+\item{\code{ptr}}{An external pointer to the underlying implementation}
+
+\item{\code{init}}{A logical variable tracking if the query condition object has been
+initialized}
+}}
+
diff --git a/apis/r/src/RcppExports.cpp b/apis/r/src/RcppExports.cpp
index 7c1121cafc..5bb38b57a2 100644
--- a/apis/r/src/RcppExports.cpp
+++ b/apis/r/src/RcppExports.cpp
@@ -553,16 +553,18 @@ BEGIN_RCPP
 END_RCPP
 }
 // resize
-void resize(const std::string& uri, Rcpp::NumericVector new_shape, std::string function_name_for_messages, Rcpp::XPtr<somactx_wrap_t> ctxxp);
-RcppExport SEXP _tiledbsoma_resize(SEXP uriSEXP, SEXP new_shapeSEXP, SEXP function_name_for_messagesSEXP, SEXP ctxxpSEXP) {
+std::string resize(const std::string& uri, Rcpp::NumericVector new_shape, std::string function_name_for_messages, bool check_only, Rcpp::XPtr<somactx_wrap_t> ctxxp);
+RcppExport SEXP _tiledbsoma_resize(SEXP uriSEXP, SEXP new_shapeSEXP, SEXP function_name_for_messagesSEXP, SEXP check_onlySEXP, SEXP ctxxpSEXP) {
 BEGIN_RCPP
+    Rcpp::RObject rcpp_result_gen;
     Rcpp::RNGScope rcpp_rngScope_gen;
     Rcpp::traits::input_parameter< const std::string& >::type uri(uriSEXP);
     Rcpp::traits::input_parameter< Rcpp::NumericVector >::type new_shape(new_shapeSEXP);
     Rcpp::traits::input_parameter< std::string >::type function_name_for_messages(function_name_for_messagesSEXP);
+    Rcpp::traits::input_parameter< bool >::type check_only(check_onlySEXP);
     Rcpp::traits::input_parameter< Rcpp::XPtr<somactx_wrap_t> >::type ctxxp(ctxxpSEXP);
-    resize(uri, new_shape, function_name_for_messages, ctxxp);
-    return R_NilValue;
+    rcpp_result_gen = Rcpp::wrap(resize(uri, new_shape, function_name_for_messages, check_only, ctxxp));
+    return rcpp_result_gen;
 END_RCPP
 }
 // resize_soma_joinid_shape
@@ -579,31 +581,35 @@ BEGIN_RCPP
 END_RCPP
 }
 // tiledbsoma_upgrade_shape
-void tiledbsoma_upgrade_shape(const std::string& uri, Rcpp::NumericVector new_shape, std::string function_name_for_messages, Rcpp::XPtr<somactx_wrap_t> ctxxp);
-RcppExport SEXP _tiledbsoma_tiledbsoma_upgrade_shape(SEXP uriSEXP, SEXP new_shapeSEXP, SEXP function_name_for_messagesSEXP, SEXP ctxxpSEXP) {
+std::string tiledbsoma_upgrade_shape(const std::string& uri, Rcpp::NumericVector new_shape, std::string function_name_for_messages, bool check_only, Rcpp::XPtr<somactx_wrap_t> ctxxp);
+RcppExport SEXP _tiledbsoma_tiledbsoma_upgrade_shape(SEXP uriSEXP, SEXP new_shapeSEXP, SEXP function_name_for_messagesSEXP, SEXP check_onlySEXP, SEXP ctxxpSEXP) {
 BEGIN_RCPP
+    Rcpp::RObject rcpp_result_gen;
     Rcpp::RNGScope rcpp_rngScope_gen;
     Rcpp::traits::input_parameter< const std::string& >::type uri(uriSEXP);
     Rcpp::traits::input_parameter< Rcpp::NumericVector >::type new_shape(new_shapeSEXP);
     Rcpp::traits::input_parameter< std::string >::type function_name_for_messages(function_name_for_messagesSEXP);
+    Rcpp::traits::input_parameter< bool >::type check_only(check_onlySEXP);
     Rcpp::traits::input_parameter< Rcpp::XPtr<somactx_wrap_t> >::type ctxxp(ctxxpSEXP);
-    tiledbsoma_upgrade_shape(uri, new_shape, function_name_for_messages, ctxxp);
-    return R_NilValue;
+    rcpp_result_gen = Rcpp::wrap(tiledbsoma_upgrade_shape(uri, new_shape, function_name_for_messages, check_only, ctxxp));
+    return rcpp_result_gen;
 END_RCPP
 }
 // upgrade_or_change_domain
-void upgrade_or_change_domain(const std::string& uri, bool is_change_domain, naxpArray nadimap, naxpSchema nadimsp, std::string function_name_for_messages, Rcpp::XPtr<somactx_wrap_t> ctxxp);
-RcppExport SEXP _tiledbsoma_upgrade_or_change_domain(SEXP uriSEXP, SEXP is_change_domainSEXP, SEXP nadimapSEXP, SEXP nadimspSEXP, SEXP function_name_for_messagesSEXP, SEXP ctxxpSEXP) {
+std::string upgrade_or_change_domain(const std::string& uri, bool is_change_domain, naxpArray nadimap, naxpSchema nadimsp, std::string function_name_for_messages, bool check_only, Rcpp::XPtr<somactx_wrap_t> ctxxp);
+RcppExport SEXP _tiledbsoma_upgrade_or_change_domain(SEXP uriSEXP, SEXP is_change_domainSEXP, SEXP nadimapSEXP, SEXP nadimspSEXP, SEXP function_name_for_messagesSEXP, SEXP check_onlySEXP, SEXP ctxxpSEXP) {
 BEGIN_RCPP
+    Rcpp::RObject rcpp_result_gen;
     Rcpp::RNGScope rcpp_rngScope_gen;
     Rcpp::traits::input_parameter< const std::string& >::type uri(uriSEXP);
     Rcpp::traits::input_parameter< bool >::type is_change_domain(is_change_domainSEXP);
     Rcpp::traits::input_parameter< naxpArray >::type nadimap(nadimapSEXP);
     Rcpp::traits::input_parameter< naxpSchema >::type nadimsp(nadimspSEXP);
     Rcpp::traits::input_parameter< std::string >::type function_name_for_messages(function_name_for_messagesSEXP);
+    Rcpp::traits::input_parameter< bool >::type check_only(check_onlySEXP);
     Rcpp::traits::input_parameter< Rcpp::XPtr<somactx_wrap_t> >::type ctxxp(ctxxpSEXP);
-    upgrade_or_change_domain(uri, is_change_domain, nadimap, nadimsp, function_name_for_messages, ctxxp);
-    return R_NilValue;
+    rcpp_result_gen = Rcpp::wrap(upgrade_or_change_domain(uri, is_change_domain, nadimap, nadimsp, function_name_for_messages, check_only, ctxxp));
+    return rcpp_result_gen;
 END_RCPP
 }
 // c_update_dataframe_schema
@@ -834,10 +840,10 @@ static const R_CallMethodDef CallEntries[] = {
     {"_tiledbsoma_c_dimnames", (DL_FUNC) &_tiledbsoma_c_dimnames, 2},
     {"_tiledbsoma_c_attrnames", (DL_FUNC) &_tiledbsoma_c_attrnames, 2},
     {"_tiledbsoma_c_schema", (DL_FUNC) &_tiledbsoma_c_schema, 2},
-    {"_tiledbsoma_resize", (DL_FUNC) &_tiledbsoma_resize, 4},
+    {"_tiledbsoma_resize", (DL_FUNC) &_tiledbsoma_resize, 5},
     {"_tiledbsoma_resize_soma_joinid_shape", (DL_FUNC) &_tiledbsoma_resize_soma_joinid_shape, 4},
-    {"_tiledbsoma_tiledbsoma_upgrade_shape", (DL_FUNC) &_tiledbsoma_tiledbsoma_upgrade_shape, 4},
-    {"_tiledbsoma_upgrade_or_change_domain", (DL_FUNC) &_tiledbsoma_upgrade_or_change_domain, 6},
+    {"_tiledbsoma_tiledbsoma_upgrade_shape", (DL_FUNC) &_tiledbsoma_tiledbsoma_upgrade_shape, 5},
+    {"_tiledbsoma_upgrade_or_change_domain", (DL_FUNC) &_tiledbsoma_upgrade_or_change_domain, 7},
     {"_tiledbsoma_c_update_dataframe_schema", (DL_FUNC) &_tiledbsoma_c_update_dataframe_schema, 6},
     {"_tiledbsoma_sr_setup", (DL_FUNC) &_tiledbsoma_sr_setup, 10},
     {"_tiledbsoma_sr_complete", (DL_FUNC) &_tiledbsoma_sr_complete, 1},
diff --git a/apis/r/src/rinterface.cpp b/apis/r/src/rinterface.cpp
index 4b23663a90..eeda0f9b22 100644
--- a/apis/r/src/rinterface.cpp
+++ b/apis/r/src/rinterface.cpp
@@ -418,18 +418,29 @@ SEXP c_schema(const std::string& uri, Rcpp::XPtr<somactx_wrap_t> ctxxp) {
 }
 
 // [[Rcpp::export]]
-void resize(
+std::string resize(
     const std::string& uri,
     Rcpp::NumericVector new_shape,
     std::string function_name_for_messages,
+    bool check_only,
     Rcpp::XPtr<somactx_wrap_t> ctxxp) {
     // This function is solely for SparseNDArray and DenseNDArray for which the
     // dims are required by the SOMA spec to be of type int64. Domain-resize for
     // variant-indexed dataframes is via upgrade_domain and change_domain.
     auto sr = tdbs::SOMAArray::open(OpenMode::write, uri, ctxxp->ctxptr);
     std::vector<int64_t> new_shape_i64 = i64_from_rcpp_numeric(new_shape);
-    sr->resize(new_shape_i64, function_name_for_messages);
+
+    std::string retval = "";
+    if (check_only) {
+        auto status_and_reason = sr->can_resize(
+            new_shape_i64, function_name_for_messages);
+        retval = status_and_reason.second;
+    } else {
+        sr->resize(new_shape_i64, function_name_for_messages);
+    }
+
     sr->close();
+    return retval;
 }
 
 // [[Rcpp::export]]
@@ -446,29 +457,40 @@ void resize_soma_joinid_shape(
 }
 
 // [[Rcpp::export]]
-void tiledbsoma_upgrade_shape(
+std::string tiledbsoma_upgrade_shape(
     const std::string& uri,
     Rcpp::NumericVector new_shape,
     std::string function_name_for_messages,
+    bool check_only,
     Rcpp::XPtr<somactx_wrap_t> ctxxp) {
     // This function is solely for SparseNDArray and DenseNDArray for which the
     // dims are required by the SOMA spec to be of type int64. Domain-resize for
     // variant-indexed dataframes is via upgrade_domain and change_domain.
     auto sr = tdbs::SOMAArray::open(OpenMode::write, uri, ctxxp->ctxptr);
     std::vector<int64_t> new_shape_i64 = i64_from_rcpp_numeric(new_shape);
-    sr->upgrade_shape(new_shape_i64, function_name_for_messages);
+
+    std::string retval = "";
+    if (check_only) {
+        auto status_and_reason = sr->can_upgrade_shape(
+            new_shape_i64, function_name_for_messages);
+        retval = status_and_reason.second;
+    } else {
+        sr->upgrade_shape(new_shape_i64, function_name_for_messages);
+    }
+
     sr->close();
+    return retval;
 }
 
 // [[Rcpp::export]]
-void upgrade_or_change_domain(
+std::string upgrade_or_change_domain(
     const std::string& uri,
     bool is_change_domain,
     naxpArray nadimap,
     naxpSchema nadimsp,
     std::string function_name_for_messages,
+    bool check_only,
     Rcpp::XPtr<somactx_wrap_t> ctxxp) {
-
     // This is pointer manipulation from R -> Rcpp SEXP -> libtiledbsoma:
     nanoarrow::UniqueArray apdim{nanoarrow_array_from_xptr(nadimap)};
     nanoarrow::UniqueSchema spdim{nanoarrow_schema_from_xptr(nadimsp)};
@@ -482,13 +504,27 @@ void upgrade_or_change_domain(
     tdbs::ArrowTable arrow_table(std::move(dimarr), std::move(dimsch));
 
     // Now call libtiledbsoma
+    std::string reason_string = "";
     auto sr = tdbs::SOMADataFrame::open(uri, OpenMode::write, ctxxp->ctxptr);
     if (is_change_domain) {
-      sr->change_domain(arrow_table, function_name_for_messages);
+        if (check_only) {
+            auto status_and_reason = sr->can_change_domain(
+                arrow_table, function_name_for_messages);
+            reason_string = status_and_reason.second;
+        } else {
+            sr->change_domain(arrow_table, function_name_for_messages);
+        }
     } else {
-      sr->upgrade_domain(arrow_table, function_name_for_messages);
+        if (check_only) {
+            auto status_and_reason = sr->can_upgrade_domain(
+                arrow_table, function_name_for_messages);
+            reason_string = status_and_reason.second;
+        } else {
+            sr->upgrade_domain(arrow_table, function_name_for_messages);
+        }
     }
     sr->close();
+    return reason_string;
 }
 
 // [[Rcpp::export]]
diff --git a/apis/r/tools/get_tarball.R b/apis/r/tools/get_tarball.R
index e08ab42d7e..7b29b5ac75 100644
--- a/apis/r/tools/get_tarball.R
+++ b/apis/r/tools/get_tarball.R
@@ -14,14 +14,14 @@ isLinux <- Sys.info()["sysname"] == "Linux"
 if (isMac) {
     arch <- system('uname -m', intern = TRUE)
     if (arch == "x86_64") {
-      url <- "https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0-rc3/tiledb-macos-x86_64-2.27.0-rc3-8d581f2.tar.gz"
+      url <- "https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0.rc5/tiledb-macos-x86_64-2.27.0.rc5-2862c30.tar.gz"
     } else if (arch == "arm64") {
-      url <- "https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0-rc3/tiledb-macos-arm64-2.27.0-rc3-8d581f2.tar.gz"
+      url <- "https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0.rc5/tiledb-macos-arm64-2.27.0.rc5-2862c30.tar.gz"
     } else {
       stop("Unsupported Mac architecture. Please have TileDB Core installed locally.")
     }
 } else if (isLinux) {
-    url <- "https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0-rc3/tiledb-linux-x86_64-2.27.0-rc3-8d581f2.tar.gz"
+    url <- "https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0.rc5/tiledb-linux-x86_64-2.27.0.rc5-2862c30.tar.gz"
 } else {
     message("Unsupported platform for downloading artifacts. Please have TileDB Core installed locally.")
     q(save = "no", status = 1)
diff --git a/apis/r/vignettes/soma-shapes.Rmd b/apis/r/vignettes/soma-shapes.Rmd
new file mode 100644
index 0000000000..bf2879fb35
--- /dev/null
+++ b/apis/r/vignettes/soma-shapes.Rmd
@@ -0,0 +1,275 @@
+---
+title: "SOMA shapes"
+output: rmarkdown::html_vignette
+vignette: >
+  %\VignetteIndexEntry{SOMA shapes}
+  %\VignetteEngine{knitr::rmarkdown}
+  %\VignetteEncoding{UTF-8}
+---
+
+```{r, include = FALSE}
+knitr::opts_chunk$set(
+  collapse = TRUE,
+  comment = "#>"
+)
+```
+
+## Overview
+
+As of TileDB-SOMA 1.15 we're proud to support a more intuitive and extensible notion of shape.
+
+Please also see the [Academy tutorial](https://cloud.tiledb.com/academy/structure/life-sciences/single-cell/tutorials/shapes/).
+
+## Example data
+
+Let's load the bundled `SOMAExperiment` containing a subsetted version of the 10X genomics [PBMC dataset](https://satijalab.github.io/seurat-object/reference/pbmc_small.html) provided by SeuratObject. This will return a `SOMAExperiment` object.
+
+```{r}
+library(tiledbsoma)
+
+exp <- load_dataset("soma-exp-pbmc-small")
+exp
+```
+
+The `obs` dataframe has a `domain`, which is a soft limit on what values can be written to it. (You'll get an error if you try to read or write `soma_joinid` values outside this range, which is an important data-integrity reassurance.)
+
+The `domain` we see here matches with the data populated inside of it.  (This will usually be the case. It might not, if you've created the dataframe but not written any data to it yet --- at that point it's empty but it still has a shape.)
+
+If you have more data --- more cells --- to add to the experiment later, you will be able resize the obs, up to the maxdomain which is a hard limit.
+
+```{r}
+exp$obs$domain()
+```
+
+```{r}
+exp$obs$maxdomain()
+```
+
+```{r}
+head(as.data.frame(exp$obs$read()$concat()))
+```
+
+The `var` dataframe's `domain` is similar:
+
+```{r}
+var <- exp$ms$get("RNA")$var
+```
+
+```{r}
+var$domain()
+```
+
+```{r}
+var$maxdomain()
+```
+
+Likewise, the N-dimensional arrays within the experiment have their shapes as well.
+
+There's an important difference: while the dataframe domain gives you the inclusive lower and upper bounds for `soma_joinid` writes, the shape for the N-dimensional arrays is the upper bound plus 1.
+
+Since there are 80 cells here and 230 genes here, `X`'s shape reflects that.
+
+```{r}
+obs <- exp$obs
+var <- exp$ms$get("RNA")$var
+X   <- exp$ms$get("RNA")$X$get("data")
+```
+
+```{r}
+obs$domain()
+```
+
+```{r}
+var$domain()
+```
+
+```{r}
+X$shape()
+```
+
+```{r}
+X$maxshape()
+```
+
+The other N-dimensional arrays are similar:
+
+```{r}
+obsm <- exp$ms$get("RNA")$obsm
+obsm$names()
+```
+
+```{r}
+obsp <- exp$ms$get("RNA")$obsp
+obsp$names()
+```
+
+```{r}
+list(
+    obsm$get("X_pca")$shape(),
+    obsm$get("X_pca")$maxshape()
+)
+```
+
+```{r}
+list(
+    obsp$get("RNA_snn")$shape(),
+    obsp$get("RNA_snn")$maxshape()
+)
+```
+
+In particular, the `X` array in this experiment --- and in most experiments --- is sparse. That means there needn't be a number in every row or cell of the matrix. Nonetheless, the shape serves as a soft limit for reads and writes: you'll get an exception trying to read or write outside of these.
+
+As a general rule you'll see the following:
+
+* An `X` array's shape is nobs x nvar
+* An `obsm` array's shape is `nobs` x some number, maybe 20
+* An `obsp` array's shape is `nobs` x `nobs`
+* A `varm` array's shape is `nvar` x some number, maybe 20
+* A `varp` array's shape is `nvar` x `nvar`
+
+## Advanced usage: dataframes with non-standard index columns
+
+In the SOMA data model, the `SOMASparseNDArray` and `SOMADenseNDArray` objects always have int64 dimensions named `soma_dim_0`, `soma_dim_1`, and up, and they have a numeric `soma_data` attribute for the contents of the array. Furthermore, this is always the case.
+
+```{r}
+exp$ms$get("RNA")$X$get("data")$schema()
+```
+
+For dataframes, though, while there must be a `soma_joinid` column of type int64, you can have one or more other index columns in addtion --- or, `soma_joinid` can be a non-index column.
+
+```{r}
+exp$obs$schema()
+```
+
+But really, dataframes are capable of more than that, via the index-column names you specify at creation time.
+
+Let's create a couple dataframes, with the same data, but different choices of index-column names.
+
+```{r}
+sdfuri1 <- withr::local_tempdir("sdf1")
+sdfuri2 <- withr::local_tempdir("sdf2")
+```
+
+```{r}
+asch <- arrow::schema(
+    arrow::field("soma_joinid", arrow::int64(), nullable = FALSE),
+    arrow::field("mystring", arrow::large_utf8(), nullable = FALSE),
+    arrow::field("myint", arrow::int32(), nullable = FALSE),
+    arrow::field("myfloat", arrow::float32(), nullable = FALSE)
+)
+
+soma_joinid = c(0, 1)
+mystring    = c("hello", "world")
+myint       = c(33, 44)
+myfloat     = c(4.5, 5.5)
+
+tbl <- arrow::arrow_table(
+    soma_joinid = c(soma_joinid),
+    mystring = c(mystring),
+    myint = c(myint),
+    myfloat = c(myfloat)
+)
+```
+
+```{r}
+sdf1 <- SOMADataFrameCreate(
+  sdfuri1,
+  asch,
+  index_column_names = c("soma_joinid", "mystring"),
+  domain = list(soma_joinid = c(0, 9), mystring = NULL)
+)
+sdf1$write(tbl)
+sdf1$close()
+```
+
+Now let's look at the `domain` and `maxdomain` for these dataframes.
+
+```{r}
+sdf1 <- SOMADataFrameOpen(sdfuri1)
+sdf1$index_column_names()
+```
+
+Here we see the `soma_joinid` slot of the dataframe's domain is as requested.
+
+Another point is that domain cannot be specified for string-type index columns.  You can set them at create one of two ways:
+
+```{r, eval=FALSE}
+domain = list(soma_joinid = (0, 9), mystring = NULL)
+```
+
+or
+
+```{r, eval=FALSE}
+domain = list(soma_joinid = (0, 9), mystring = c('', ''))
+```
+
+and in either case the domain slot for a string-typed index column will read back as `('', '')`.
+
+```{r}
+sdf1$domain()
+```
+
+```{r}
+sdf1$maxdomain()
+```
+
+Now let's look at our other dataframe. Here `soma_joinid` is not an index column at all. This is fine, as long as within the data you write to it, the index-column values uniquely identify each row.
+
+```{r}
+sdf2 <- SOMADataFrameCreate(
+  sdfuri2,
+  asch,
+  index_column_names = c("myfloat", "myint"),
+  domain = list(myfloat = c(0, 9999), myint = c(-1000, 1000))
+)
+sdf2$write(tbl)
+sdf2$close()
+```
+
+```{r}
+sdf2 <- SOMADataFrameOpen(sdfuri1)
+sdf2$index_column_names()
+```
+
+The domain reads back as written:
+
+```{r}
+sdf2$domain()
+```
+
+```{r}
+sdf2$maxdomain()
+```
+
+## Advanced usage: using resize at the dataframe/array level
+
+In the TileDB-SOMA Python API, there is a method for resizing all the dataframes and arrays within an experiment. At present we do not yet offer a corresponding method in the TileDB-SOMA R API, for the simple reason that there is low demand for it. Nonetheless, for completeness, we offer here guidance on how to resizes dataframes and arrays within a TileDB-SOMA experiment.
+
+For N-dimensional arrays that have been upgraded, or that were created using TileDB-SOMA 1.15 or higher, simply do the following:
+
+* If the array's `$tiledbsoma_has_upgraded_shape()` reports `FALSE`, invoke the `$tiledbsoma_upgrade_shape()` method.
+* Otherwise invoke the `$.resize()` method.
+
+Let's do a fresh unpack of a pre-1.15 experiment:
+
+```{r}
+exp <- load_dataset("soma-exp-pbmc-small-pre-1.15")
+exp
+```
+
+Here we see that the X array has not been upgraded, and that its shape reports the same as maxshape:
+
+```{r}
+X <- exp$ms$get("RNA")$X$get("data")
+X$tiledbsoma_has_upgraded_shape()
+```
+
+```{r}
+X$shape()
+```
+
+```{r}
+X$maxshape()
+```
+
+Given that pre-1.15 TileDB-SOMA-R arrays were created with a maxshape leaving no room for growth, these arrays cannot have their shape resized any further. From 1.15 onward, of course, as we've see above, arrays are created with room for growth and you can resize them upward.
diff --git a/doc/source/notebooks/tutorial_soma_shape.ipynb b/doc/source/notebooks/tutorial_soma_shape.ipynb
new file mode 100644
index 0000000000..eab6b9470f
--- /dev/null
+++ b/doc/source/notebooks/tutorial_soma_shape.ipynb
@@ -0,0 +1,1 @@
+../../../apis/python/notebooks/tutorial_soma_shape.ipynb
\ No newline at end of file
diff --git a/doc/source/notebooks/tutorial_spatial.ipynb b/doc/source/notebooks/tutorial_spatial.ipynb
new file mode 100644
index 0000000000..62565ac954
--- /dev/null
+++ b/doc/source/notebooks/tutorial_spatial.ipynb
@@ -0,0 +1,1 @@
+../../../apis/python/notebooks/tutorial_spatial.ipynb
\ No newline at end of file
diff --git a/libtiledbsoma/cmake/Modules/FindTileDB_EP.cmake b/libtiledbsoma/cmake/Modules/FindTileDB_EP.cmake
index 8fec48e655..07f8515f64 100644
--- a/libtiledbsoma/cmake/Modules/FindTileDB_EP.cmake
+++ b/libtiledbsoma/cmake/Modules/FindTileDB_EP.cmake
@@ -58,8 +58,8 @@ else()
     # NB When updating the pinned URLs here, please also update in file apis/r/tools/get_tarball.R
     if(DOWNLOAD_TILEDB_PREBUILT)
         if (WIN32) # Windows
-          SET(DOWNLOAD_URL "https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0-rc3/tiledb-windows-x86_64-2.27.0-rc3-8d581f2.zip")
-          SET(DOWNLOAD_SHA1 "db13fe37a28f25962655e70854ca50435121287d")
+          SET(DOWNLOAD_URL "https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0.rc5/tiledb-windows-x86_64-2.27.0.rc5-2862c30.zip")
+          SET(DOWNLOAD_SHA1 "c23d84a55dac572c767d461e7b2a697edce88512")
         elseif(APPLE) # OSX
 
           # Status quo as of 2023-05-18:
@@ -76,22 +76,22 @@ else()
           #   o CMAKE_SYSTEM_PROCESSOR is x86_64
 
           if (CMAKE_OSX_ARCHITECTURES STREQUAL x86_64)
-            SET(DOWNLOAD_URL "https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0-rc3/tiledb-macos-x86_64-2.27.0-rc3-8d581f2.tar.gz")
-            SET(DOWNLOAD_SHA1 "470c236763f731c018badf7033e0bdd376d98db0")
+            SET(DOWNLOAD_URL "https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0.rc5/tiledb-macos-x86_64-2.27.0.rc5-2862c30.tar.gz")
+            SET(DOWNLOAD_SHA1 "586020a2ee2c93eda3013cad17b4dd20e6b0b8be")
           elseif (CMAKE_OSX_ARCHITECTURES STREQUAL arm64)
-            SET(DOWNLOAD_URL "https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0-rc3/tiledb-macos-arm64-2.27.0-rc3-8d581f2.tar.gz")
-            SET(DOWNLOAD_SHA1 "4f6929cb9e1631989b0e45fb5933022bd387db31")
+            SET(DOWNLOAD_URL "https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0.rc5/tiledb-macos-arm64-2.27.0.rc5-2862c30.tar.gz")
+            SET(DOWNLOAD_SHA1 "ebcf68fbd5ba6369123d5cc4adaa37aec9440768")
           elseif (CMAKE_SYSTEM_PROCESSOR MATCHES "(x86_64)|(AMD64|amd64)|(^i.86$)")
-            SET(DOWNLOAD_URL "https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0-rc3/tiledb-macos-x86_64-2.27.0-rc3-8d581f2.tar.gz")
-            SET(DOWNLOAD_SHA1 "470c236763f731c018badf7033e0bdd376d98db0")
+            SET(DOWNLOAD_URL "https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0.rc5/tiledb-macos-x86_64-2.27.0.rc5-2862c30.tar.gz")
+            SET(DOWNLOAD_SHA1 "586020a2ee2c93eda3013cad17b4dd20e6b0b8be")
           elseif (CMAKE_SYSTEM_PROCESSOR MATCHES "^aarch64" OR CMAKE_SYSTEM_PROCESSOR MATCHES "^arm")
-            SET(DOWNLOAD_URL "https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0-rc3/tiledb-macos-arm64-2.27.0-rc3-8d581f2.tar.gz")
-            SET(DOWNLOAD_SHA1 "4f6929cb9e1631989b0e45fb5933022bd387db31")
+            SET(DOWNLOAD_URL "https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0.rc5/tiledb-macos-arm64-2.27.0.rc5-2862c30.tar.gz")
+            SET(DOWNLOAD_SHA1 "ebcf68fbd5ba6369123d5cc4adaa37aec9440768")
           endif()
 
         else() # Linux
-          SET(DOWNLOAD_URL "https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0-rc3/tiledb-linux-x86_64-2.27.0-rc3-8d581f2.tar.gz")
-          SET(DOWNLOAD_SHA1 "750e5629fcab223de947f395a8148fe0ba163af8")
+          SET(DOWNLOAD_URL "https://github.com/TileDB-Inc/TileDB/releases/download/2.27.0.rc5/tiledb-linux-x86_64-2.27.0.rc5-2862c30.tar.gz")
+          SET(DOWNLOAD_SHA1 "9c1b1810cedd18bf30d9470d4569c9273e7f0dc3")
         endif()
 
         ExternalProject_Add(ep_tiledb
@@ -113,8 +113,8 @@ else()
     else() # Build from source
         ExternalProject_Add(ep_tiledb
           PREFIX "externals"
-          URL "https://github.com/TileDB-Inc/TileDB/archive/2.27.0-rc3.zip"
-          URL_HASH SHA1=ecfb3121208715b3099e101af0b1460c5cda7510
+          URL "https://github.com/TileDB-Inc/TileDB/archive/2.27.0.rc5.zip"
+          URL_HASH SHA1=a7fae38918f268e6efbcc688d51800b93cfe3f5d
           DOWNLOAD_NAME "tiledb.zip"
           CMAKE_ARGS
             -DCMAKE_INSTALL_PREFIX=${EP_INSTALL_PREFIX}
diff --git a/libtiledbsoma/src/soma/managed_query.cc b/libtiledbsoma/src/soma/managed_query.cc
index 3345b32a5b..7aad1e67fc 100644
--- a/libtiledbsoma/src/soma/managed_query.cc
+++ b/libtiledbsoma/src/soma/managed_query.cc
@@ -202,6 +202,104 @@ void ManagedQuery::submit_read() {
     });
 }
 
+std::optional<std::shared_ptr<ArrayBuffers>> ManagedQuery::read_next() {
+    setup_read();
+
+    if (is_empty_query() && !query_submitted_) {
+        query_submitted_ = true;
+        return buffers_;
+    }
+
+    if (is_complete(false)) {
+        return std::nullopt;
+    }
+
+    query_submitted_ = true;
+    query_future_ = std::async(std::launch::async, [&]() {
+        LOG_DEBUG("[ManagedQuery] submit thread start");
+        try {
+            query_->submit();
+        } catch (const std::exception& e) {
+            return StatusAndException(false, e.what());
+        }
+        LOG_DEBUG("[ManagedQuery] submit thread done");
+        return StatusAndException(true, "success");
+    });
+
+    if (query_future_.valid()) {
+        LOG_DEBUG(fmt::format("[ManagedQuery] [{}] Waiting for query", name_));
+        query_future_.wait();
+        LOG_DEBUG(
+            fmt::format("[ManagedQuery] [{}] Done waiting for query", name_));
+
+        auto retval = query_future_.get();
+        if (!retval.succeeded()) {
+            throw TileDBSOMAError(fmt::format(
+                "[ManagedQuery] [{}] Query FAILED: {}",
+                name_,
+                retval.message()));
+        }
+
+    } else {
+        throw TileDBSOMAError(
+            fmt::format("[ManagedQuery] [{}] 'query_future_' invalid", name_));
+    }
+
+    auto status = query_->query_status();
+
+    if (status == Query::Status::FAILED) {
+        throw TileDBSOMAError(
+            fmt::format("[ManagedQuery] [{}] Query FAILED", name_));
+    }
+
+    // If the query was ever incomplete, the result buffers contents are not
+    // complete.
+    if (status == Query::Status::INCOMPLETE) {
+        results_complete_ = false;
+    } else if (status == Query::Status::COMPLETE) {
+        results_complete_ = true;
+    }
+
+    // Update ColumnBuffer size to match query results
+    size_t num_cells = 0;
+    for (auto& name : buffers_->names()) {
+        num_cells = buffers_->at(name)->update_size(*query_);
+        LOG_DEBUG(fmt::format(
+            "[ManagedQuery] [{}] Buffer {} cells={}", name_, name, num_cells));
+    }
+    total_num_cells_ += num_cells;
+
+    // TODO: retry the query with larger buffers
+    if (status == Query::Status::INCOMPLETE && !num_cells) {
+        throw TileDBSOMAError(
+            fmt::format("[ManagedQuery] [{}] Buffers are too small.", name_));
+    }
+
+    // Visit all attributes and retrieve enumeration vectors
+    auto attribute_map = schema_->attributes();
+    for (auto& nmit : attribute_map) {
+        auto attrname = nmit.first;
+        auto attribute = nmit.second;
+        auto enumname = AttributeExperimental::get_enumeration_name(
+            *ctx_, attribute);
+        if (enumname != std::nullopt) {
+            auto enumeration = ArrayExperimental::get_enumeration(
+                *ctx_, *array_, enumname.value());
+            auto enumvec = enumeration.as_vector<std::string>();
+            if (!buffers_->contains(attrname)) {
+                continue;
+            }
+            auto colbuf = buffers_->at(attrname);
+            colbuf->add_enumeration(enumvec);
+            LOG_DEBUG(fmt::format(
+                "[ManagedQuery] got Enumeration '{}' for attribute '{}'",
+                enumname.value(),
+                attrname));
+        }
+    }
+    return buffers_;
+}
+
 // Please see the header-file comments for context.
 void ManagedQuery::_fill_in_subarrays_if_dense(bool is_read) {
     LOG_TRACE("[ManagedQuery] _fill_in_subarrays enter");
@@ -380,6 +478,7 @@ void ManagedQuery::_fill_in_subarrays_if_dense_with_new_shape(
 
 std::shared_ptr<ArrayBuffers> ManagedQuery::results() {
     if (is_empty_query()) {
+        query_submitted_ = true;
         return buffers_;
     }
 
diff --git a/libtiledbsoma/src/soma/managed_query.h b/libtiledbsoma/src/soma/managed_query.h
index 2f806b0eee..ab76a48482 100644
--- a/libtiledbsoma/src/soma/managed_query.h
+++ b/libtiledbsoma/src/soma/managed_query.h
@@ -281,6 +281,8 @@ class ManagedQuery {
      */
     void setup_read();
 
+    std::optional<std::shared_ptr<ArrayBuffers>> read_next();
+
     /**
      * @brief Check if the query is complete.
      *
@@ -423,6 +425,10 @@ class ManagedQuery {
         return query_->query_status();
     }
 
+    bool is_first_read() const {
+        return !query_submitted_;
+    }
+
    private:
     //===================================================================
     //= private non-static
diff --git a/libtiledbsoma/src/soma/soma_array.h b/libtiledbsoma/src/soma/soma_array.h
index c4792ecd4a..06d5b46cb8 100644
--- a/libtiledbsoma/src/soma/soma_array.h
+++ b/libtiledbsoma/src/soma/soma_array.h
@@ -380,7 +380,7 @@ class SOMAArray : public SOMAObject {
         if (partition_index >= partition_count) {
             // TODO this use to be formatted with fmt::format which is part of
             // internal header spd/log/fmt/fmt.h and should not be used.
-            // In C++20, this can be replaced with std::format.
+            // In C++20, this can be replaced with fmt::format.
             std::ostringstream err;
             err << "[SOMAArray] partition_index (" << partition_index
                 << ") must be < partition_count (" << partition_count;
@@ -398,7 +398,7 @@ class SOMAArray : public SOMAObject {
 
             // TODO this use to be formatted with fmt::format which is part of
             // internal header spd/log/fmt/fmt.h and should not be used.
-            // In C++20, this can be replaced with std::format.
+            // In C++20, this can be replaced with fmt::format.
             std::ostringstream log_dbg;
             log_dbg << "[SOMAArray] set_dim_points partitioning:"
                     << " sizeof(T)=" << sizeof(T) << " dim=" << dim
diff --git a/libtiledbsoma/src/soma/soma_scene.cc b/libtiledbsoma/src/soma/soma_scene.cc
index a0cd5b77ac..680d21a186 100644
--- a/libtiledbsoma/src/soma/soma_scene.cc
+++ b/libtiledbsoma/src/soma/soma_scene.cc
@@ -46,30 +46,10 @@ void SOMAScene::create(
     std::optional<TimestampRange> timestamp) {
     try {
         std::filesystem::path scene_uri(uri);
-
-        SOMAGroup::create(ctx, scene_uri.string(), "SOMAScene", timestamp);
-        SOMACollection::create((scene_uri / "img").string(), ctx, timestamp);
-        SOMACollection::create((scene_uri / "obsl").string(), ctx, timestamp);
-        SOMACollection::create((scene_uri / "varl").string(), ctx, timestamp);
-
-        auto name = std::string(std::filesystem::path(uri).filename());
-        auto group = SOMAGroup::open(
-            OpenMode::write, scene_uri.string(), ctx, name, timestamp);
-        group->set(
-            (scene_uri / "img").string(),
-            URIType::absolute,
-            "img",
-            "SOMACollection");
-        group->set(
-            (scene_uri / "obsl").string(),
-            URIType::absolute,
-            "obsl",
-            "SOMACollection");
-        group->set(
-            (scene_uri / "varl").string(),
-            URIType::absolute,
-            "varl",
-            "SOMAColelction");
+        auto group = SOMAGroup::create(
+            ctx, scene_uri.string(), "SOMAScene", timestamp);
+        // TODO: Set extra metadata.
+        // (See story: https://app.shortcut.com/tiledb-inc/story/59915)
         group->close();
     } catch (TileDBError& e) {
         throw TileDBSOMAError(e.what());
diff --git a/libtiledbsoma/src/soma/soma_scene.h b/libtiledbsoma/src/soma/soma_scene.h
index 6ac36df3c7..f056961f65 100644
--- a/libtiledbsoma/src/soma/soma_scene.h
+++ b/libtiledbsoma/src/soma/soma_scene.h
@@ -51,7 +51,7 @@ class SOMAScene : public SOMACollection {
      *
      * @param uri URI to create the SOMAScene
      * @param schema TileDB ArraySchema
-     * @param platform_config Optional config parameter dictionary
+     * @param timestamp Optional pair indicating timestamp start and end
      */
     static void create(
         std::string_view uri,


Test patch:
diff --git a/apis/python/tests/test_basic_anndata_io.py b/apis/python/tests/test_basic_anndata_io.py
index 505394bb06..7805cbe0f6 100644
--- a/apis/python/tests/test_basic_anndata_io.py
+++ b/apis/python/tests/test_basic_anndata_io.py
@@ -7,6 +7,7 @@
 import anndata
 import numpy as np
 import pandas as pd
+import pyarrow as pa
 import pytest
 import scipy
 import somacore
@@ -1347,3 +1348,134 @@ def test_nan_append(conftest_pbmc_small, dtype, nans, new_obs_ids):
         measurement_name="RNA",
         registration_mapping=rd,
     )
+
+
+# This tests https://github.com/single-cell-data/TileDB-SOMA/pull/3354:
+# * High-cardinality categorical-of-string gets converted to plain string, for
+#   a significant performance benefit in subsequent processing
+# * Two particular case needing testing:
+#   o When a high-cardinality enum is appended to existing storage with categorical of string;
+#   o When a low-cardinality enum is appended to existing storage with plain string.
+def test_decat_append(tmp_path):
+
+    # Prepare the AnnData inputs
+    nobs_under = tiledbsoma.io.conversions.STRING_DECAT_THRESHOLD - 2
+    nobs_over = tiledbsoma.io.conversions.STRING_DECAT_THRESHOLD + 2
+    nvar = 100
+
+    obs_ids_under = [f"under_{e:08}" for e in range(nobs_under)]
+    obs_ids_over = [f"over_{e:08}" for e in range(nobs_over)]
+    var_ids = [f"gene_{e:08}" for e in range(nvar)]
+
+    enum_values_under = [f"enum_u_{e:06}" for e in range(nobs_under)]
+    enum_values_over = [f"enum_o_{e:06}" for e in range(nobs_over)]
+
+    obs_under = pd.DataFrame(
+        data={
+            "obs_id": np.asarray(obs_ids_under),
+            "is_primary_data": np.asarray([True] * nobs_under),
+            "myenum": pd.Series(np.asarray(enum_values_under), dtype="category"),
+        }
+    )
+    obs_under.set_index("obs_id", inplace=True)
+
+    obs_over = pd.DataFrame(
+        data={
+            "obs_id": np.asarray(obs_ids_over),
+            "is_primary_data": np.asarray([True] * nobs_over),
+            "myenum": pd.Series(np.asarray(enum_values_over), dtype="category"),
+        }
+    )
+    obs_over.set_index("obs_id", inplace=True)
+
+    var = pd.DataFrame(
+        data={
+            "var_id": np.asarray(var_ids),
+            "mybool": np.asarray([True] * nvar),
+        }
+    )
+    var.set_index("var_id", inplace=True)
+
+    X_under = scipy.sparse.random(
+        nobs_under, nvar, density=0.1, dtype=np.float64
+    ).tocsr()
+    X_over = scipy.sparse.random(nobs_over, nvar, density=0.1, dtype=np.float64).tocsr()
+
+    adata_under = anndata.AnnData(
+        X=X_under, obs=obs_under, var=var, dtype=X_under.dtype
+    )
+    adata_over = anndata.AnnData(X=X_over, obs=obs_over, var=var, dtype=X_over.dtype)
+
+    # Do the initial ingests from AnnData format to TileDB-SOMA format
+    path_under = (tmp_path / "under").as_posix()
+    path_over = (tmp_path / "over").as_posix()
+
+    tiledbsoma.io.from_anndata(path_under, adata_under, "RNA")
+    tiledbsoma.io.from_anndata(path_over, adata_over, "RNA")
+
+    # Check that the low-cardinality categorical-of-string in the AnnData has
+    # been ingested to TileDB-SOMA enum-of-string.
+    with tiledbsoma.Experiment.open(path_under) as exp_under:
+        assert pa.types.is_dictionary(exp_under.obs.schema.field("myenum").type)
+        obs_table = exp_under.obs.read().concat()
+        assert obs_table.column("myenum").to_pylist() == enum_values_under
+
+    # Check that the high-cardinality categorical-of-string in the AnnData has
+    # been ingested to TileDB-SOMA plain string.
+    with tiledbsoma.Experiment.open(path_over) as exp_over:
+        assert not pa.types.is_dictionary(exp_over.obs.schema.field("myenum").type)
+        obs_table = exp_over.obs.read().concat()
+        assert obs_table.column("myenum").to_pylist() == enum_values_over
+
+    # Append over-the-threshold AnnData to under-the-threshold TileDB-SOMA
+    # storage, and vice versa.
+    rd_under_over = tiledbsoma.io.register_anndatas(
+        experiment_uri=path_under,
+        adatas=[adata_over],
+        measurement_name="RNA",
+        obs_field_name="obs_id",
+        var_field_name="var_id",
+    )
+
+    rd_over_under = tiledbsoma.io.register_anndatas(
+        experiment_uri=path_over,
+        adatas=[adata_under],
+        measurement_name="RNA",
+        obs_field_name="obs_id",
+        var_field_name="var_id",
+    )
+
+    tiledbsoma.io.resize_experiment(
+        path_under,
+        nobs=rd_under_over.get_obs_shape(),
+        nvars=rd_under_over.get_var_shapes(),
+    )
+    tiledbsoma.io.resize_experiment(
+        path_over,
+        nobs=rd_over_under.get_obs_shape(),
+        nvars=rd_over_under.get_var_shapes(),
+    )
+
+    tiledbsoma.io.from_anndata(
+        path_under, adata_over, "RNA", registration_mapping=rd_under_over
+    )
+    tiledbsoma.io.from_anndata(
+        path_over, adata_under, "RNA", registration_mapping=rd_over_under
+    )
+
+    # Check that the appends happened successfully
+    with tiledbsoma.Experiment.open(path_under) as exp_under:
+        assert pa.types.is_dictionary(exp_under.obs.schema.field("myenum").type)
+        obs_table = exp_under.obs.read().concat()
+        assert (
+            obs_table.column("myenum").to_pylist()
+            == enum_values_under + enum_values_over
+        )
+
+    with tiledbsoma.Experiment.open(path_over) as exp_over:
+        assert not pa.types.is_dictionary(exp_over.obs.schema.field("myenum").type)
+        obs_table = exp_over.obs.read().concat()
+        assert (
+            obs_table.column("myenum").to_pylist()
+            == enum_values_over + enum_values_under
+        )
diff --git a/apis/python/tests/test_basic_xarray_io.py b/apis/python/tests/test_basic_xarray_io.py
index c04b04d51b..8bac6df84f 100644
--- a/apis/python/tests/test_basic_xarray_io.py
+++ b/apis/python/tests/test_basic_xarray_io.py
@@ -71,15 +71,16 @@ def test_getitem(self, xr_soma_data_array, xr_numpy_data_array, key):
 
         np.testing.assert_equal(actual.data.compute(), expected.data)
 
-    @pytest.mark.parametrize(
-        "key",
-        [
-            (slice(0, 4, 2),),
-        ],
-    )
-    def test_getitem_with_steps(self, xr_soma_data_array, key):
-        with pytest.raises(NotImplementedError):
-            xr_soma_data_array[key].data.compute()
+    # TEMPORARY PENDING https://github.com/single-cell-data/TileDB-SOMA/issues/3398
+    # @pytest.mark.parametrize(
+    #    "key",
+    #    [
+    #        (slice(0, 4, 2),),
+    #    ],
+    # )
+    # def test_getitem_with_steps(self, xr_soma_data_array, key):
+    #    with pytest.raises(NotImplementedError):
+    #        xr_soma_data_array[key].data.compute()
 
 
 class TestDenseNDDataArray3D:
@@ -151,12 +152,13 @@ def test_getitem(self, xr_soma_data_array, xr_numpy_data_array, key):
 
         np.testing.assert_equal(actual.data.compute(), expected.data)
 
-    @pytest.mark.parametrize(
-        "key",
-        [
-            (0, 1, slice(0, 4, 2)),
-        ],
-    )
-    def test_getitem_with_steps(self, xr_soma_data_array, key):
-        with pytest.raises(NotImplementedError):
-            xr_soma_data_array[key].data.compute()
+    # TEMPORARY PENDING https://github.com/single-cell-data/TileDB-SOMA/issues/3398
+    # @pytest.mark.parametrize(
+    #    "key",
+    #    [
+    #        (0, 1, slice(0, 4, 2)),
+    #    ],
+    # )
+    # def test_getitem_with_steps(self, xr_soma_data_array, key):
+    #    with pytest.raises(NotImplementedError):
+    #        xr_soma_data_array[key].data.compute()
diff --git a/apis/python/tests/test_context.py b/apis/python/tests/test_context.py
index 14d474abe5..22a20e2004 100644
--- a/apis/python/tests/test_context.py
+++ b/apis/python/tests/test_context.py
@@ -1,4 +1,5 @@
 import datetime
+import importlib
 import time
 from unittest import mock
 
@@ -122,3 +123,40 @@ def test_malformed_concurrency_config_value():
         tiledbsoma.IntIndexer(np.arange(100, dtype=np.int64), context=ctx).get_indexer(
             np.array([0, 1])
         )
+
+
+@pytest.mark.skipif(
+    importlib.util.find_spec("tiledb") is not None, reason="TileDB-Py is installed"
+)
+def test_tiledb_ctx_without_tiledb():
+    # Test that tiledb_ctx errors out as expected without tiledb-py
+
+    with pytest.raises(ModuleNotFoundError):
+        tiledbsoma.SOMATileDBContext(tiledb_ctx="junk")
+
+    sctx = tiledbsoma.SOMATileDBContext()
+    with pytest.raises(ModuleNotFoundError):
+        sctx.tiledb_ctx
+
+    with pytest.raises(ModuleNotFoundError):
+        sctx.replace(tiledb_ctx="junk")
+
+
+@pytest.mark.skipif(
+    importlib.util.find_spec("tiledb") is None, reason="TileDB-Py is not installed"
+)
+def test_tiledb_ctx_with_tiledb():
+    # If tiledb-py is installed, test that tiledb_ctx works to handle tiledb.Ctx
+    import tiledb
+
+    # Default
+    sctx = tiledbsoma.SOMATileDBContext(tiledb_ctx=tiledb.Ctx())
+    assert sctx.tiledb_ctx.config() == tiledb.Ctx().config()
+
+    # Pass config
+    sctx = tiledbsoma.SOMATileDBContext(tiledb_ctx=tiledb.Ctx({"foo": "bar"}))
+    assert sctx.tiledb_ctx.config()["foo"] == "bar"
+
+    # Replace config
+    sctx = sctx.replace(tiledb_ctx=tiledb.Ctx({"foo": "baz"}))
+    assert sctx.tiledb_ctx.config()["foo"] == "baz"
diff --git a/apis/python/tests/test_dataframe.py b/apis/python/tests/test_dataframe.py
index df2dff0a50..c607cfcf7f 100644
--- a/apis/python/tests/test_dataframe.py
+++ b/apis/python/tests/test_dataframe.py
@@ -21,10 +21,10 @@ def arrow_schema():
     def _schema():
         return pa.schema(
             [
-                pa.field("foo", pa.int64()),
-                pa.field("bar", pa.float64()),
-                pa.field("baz", pa.string()),
-                pa.field("quux", pa.bool_()),
+                pa.field("myint", pa.int64()),
+                pa.field("myfloat", pa.float64()),
+                pa.field("mystring", pa.string()),
+                pa.field("mybool", pa.bool_()),
             ]
         )
 
@@ -36,10 +36,10 @@ def test_dataframe(tmp_path, arrow_schema):
 
     asch = pa.schema(
         [
-            ("foo", pa.int32()),
-            ("bar", pa.float64()),
-            ("baz", pa.large_string()),
-            ("quux", pa.bool_()),
+            ("myint", pa.int32()),
+            ("myfloat", pa.float64()),
+            ("mystring", pa.large_string()),
+            ("mybool", pa.bool_()),
         ]
     )
 
@@ -55,7 +55,7 @@ def test_dataframe(tmp_path, arrow_schema):
         # nonexistent indexed column
         soma.DataFrame.create(uri, schema=asch, index_column_names=["bogus"])
     soma.DataFrame.create(
-        uri, schema=asch, index_column_names=["foo"], domain=[[0, 99]]
+        uri, schema=asch, index_column_names=["myint"], domain=[[0, 99]]
     ).close()
 
     assert soma.DataFrame.exists(uri)
@@ -67,7 +67,7 @@ def test_dataframe(tmp_path, arrow_schema):
         assert len(sdf) == 0
 
         assert sorted(sdf.schema.names) == sorted(
-            ["foo", "bar", "baz", "soma_joinid", "quux"]
+            ["myint", "myfloat", "mystring", "soma_joinid", "mybool"]
         )
         assert sorted(sdf.keys()) == sorted(sdf.schema.names)
 
@@ -76,10 +76,10 @@ def test_dataframe(tmp_path, arrow_schema):
         for _ in range(3):
             pydict = {}
             pydict["soma_joinid"] = [0, 1, 2, 3, 4]
-            pydict["foo"] = [10, 20, 30, 40, 50]
-            pydict["bar"] = [4.1, 5.2, 6.3, 7.4, 8.5]
-            pydict["baz"] = ["apple", "ball", "cat", "dog", "egg"]
-            pydict["quux"] = [True, False, False, True, False]
+            pydict["myint"] = [10, 20, 30, 40, 50]
+            pydict["myfloat"] = [4.1, 5.2, 6.3, 7.4, 8.5]
+            pydict["mystring"] = ["apple", "ball", "cat", "dog", "egg"]
+            pydict["mybool"] = [True, False, False, True, False]
             rb = pa.Table.from_pydict(pydict)
 
             sdf.tiledbsoma_resize_soma_joinid_shape(len(rb))
@@ -108,20 +108,20 @@ def test_dataframe(tmp_path, arrow_schema):
         assert table.num_rows == 5
         assert table.num_columns == 5
         assert [e.as_py() for e in table["soma_joinid"]] == pydict["soma_joinid"]
-        assert [e.as_py() for e in table["foo"]] == pydict["foo"]
-        assert [e.as_py() for e in table["bar"]] == pydict["bar"]
-        assert [e.as_py() for e in table["baz"]] == pydict["baz"]
-        assert [e.as_py() for e in table["quux"]] == pydict["quux"]
+        assert [e.as_py() for e in table["myint"]] == pydict["myint"]
+        assert [e.as_py() for e in table["myfloat"]] == pydict["myfloat"]
+        assert [e.as_py() for e in table["mystring"]] == pydict["mystring"]
+        assert [e.as_py() for e in table["mybool"]] == pydict["mybool"]
 
         # Read ids
         table = sdf.read(coords=[[30, 10]]).concat()
         assert table.num_rows == 2
         assert table.num_columns == 5
         assert sorted([e.as_py() for e in table["soma_joinid"]]) == [0, 2]
-        assert sorted([e.as_py() for e in table["foo"]]) == [10, 30]
-        assert sorted([e.as_py() for e in table["bar"]]) == [4.1, 6.3]
-        assert sorted([e.as_py() for e in table["baz"]]) == ["apple", "cat"]
-        assert [e.as_py() for e in table["quux"]] == [True, False]
+        assert sorted([e.as_py() for e in table["myint"]]) == [10, 30]
+        assert sorted([e.as_py() for e in table["myfloat"]]) == [4.1, 6.3]
+        assert sorted([e.as_py() for e in table["mystring"]]) == ["apple", "cat"]
+        assert [e.as_py() for e in table["mybool"]] == [True, False]
 
     # Open and read with bindings
     with contextlib.closing(
@@ -133,18 +133,18 @@ def test_dataframe(tmp_path, arrow_schema):
         assert table.num_rows == 5
         assert table.num_columns == 5
         assert [e.as_py() for e in table["soma_joinid"]] == pydict["soma_joinid"]
-        assert [e.as_py() for e in table["foo"]] == pydict["foo"]
-        assert [e.as_py() for e in table["bar"]] == pydict["bar"]
-        assert [e.as_py() for e in table["baz"]] == pydict["baz"]
-        assert [e.as_py() for e in table["quux"]] == pydict["quux"]
+        assert [e.as_py() for e in table["myint"]] == pydict["myint"]
+        assert [e.as_py() for e in table["myfloat"]] == pydict["myfloat"]
+        assert [e.as_py() for e in table["mystring"]] == pydict["mystring"]
+        assert [e.as_py() for e in table["mybool"]] == pydict["mybool"]
 
     with soma.DataFrame.open(uri) as A:
         cfg = A.config_options_from_schema()
         assert not cfg.allows_duplicates
-        assert json.loads(cfg.dims)["foo"]["filters"] == [
+        assert json.loads(cfg.dims)["myint"]["filters"] == [
             {"COMPRESSION_LEVEL": 3, "name": "ZSTD"}
         ]
-        assert json.loads(cfg.attrs)["bar"]["filters"] == [
+        assert json.loads(cfg.attrs)["myfloat"]["filters"] == [
             {"COMPRESSION_LEVEL": -1, "name": "ZSTD"}
         ]
 
@@ -202,16 +202,16 @@ def test_dataframe_reopen(tmp_path, arrow_schema):
 
 def test_dataframe_with_float_dim(tmp_path, arrow_schema):
     sdf = soma.DataFrame.create(
-        tmp_path.as_posix(), schema=arrow_schema(), index_column_names=("bar",)
+        tmp_path.as_posix(), schema=arrow_schema(), index_column_names=("myfloat",)
     )
-    assert sdf.index_column_names == ("bar",)
+    assert sdf.index_column_names == ("myfloat",)
 
 
 def test_dataframe_with_enumeration(tmp_path):
     schema = pa.schema(
         [
-            pa.field("foo", pa.dictionary(pa.int64(), pa.large_string())),
-            pa.field("bar", pa.dictionary(pa.int64(), pa.large_string())),
+            pa.field("myint", pa.dictionary(pa.int64(), pa.large_string())),
+            pa.field("myfloat", pa.dictionary(pa.int64(), pa.large_string())),
         ]
     )
     enums = {"enmr1": ("a", "bb", "ccc"), "enmr2": ("cat", "dog")}
@@ -220,19 +220,19 @@ def test_dataframe_with_enumeration(tmp_path):
     ) as sdf:
         data = {}
         data["soma_joinid"] = [0, 1, 2, 3, 4]
-        data["foo"] = ["a", "bb", "ccc", "bb", "a"]
-        data["bar"] = ["cat", "dog", "cat", "cat", "cat"]
+        data["myint"] = ["a", "bb", "ccc", "bb", "a"]
+        data["myfloat"] = ["cat", "dog", "cat", "cat", "cat"]
         with pytest.raises(soma.SOMAError):
             sdf.write(pa.Table.from_pydict(data))
 
-        data["foo"] = pd.Categorical(["a", "bb", "ccc", "bb", "a"])
-        data["bar"] = pd.Categorical(["cat", "dog", "cat", "cat", "cat"])
+        data["myint"] = pd.Categorical(["a", "bb", "ccc", "bb", "a"])
+        data["myfloat"] = pd.Categorical(["cat", "dog", "cat", "cat", "cat"])
         sdf.write(pa.Table.from_pydict(data))
 
     with soma.DataFrame.open(tmp_path.as_posix()) as sdf:
         df = sdf.read().concat()
-        np.testing.assert_array_equal(df["foo"].chunk(0).dictionary, enums["enmr1"])
-        np.testing.assert_array_equal(df["bar"].chunk(0).dictionary, enums["enmr2"])
+        np.testing.assert_array_equal(df["myint"].chunk(0).dictionary, enums["enmr1"])
+        np.testing.assert_array_equal(df["myfloat"].chunk(0).dictionary, enums["enmr2"])
 
 
 @pytest.fixture
@@ -1749,8 +1749,8 @@ def test_only_evolve_schema_when_enmr_is_extended(tmp_path):
 
     schema = pa.schema(
         [
-            pa.field("foo", pa.dictionary(pa.int64(), pa.large_string())),
-            pa.field("bar", pa.large_string()),
+            pa.field("myint", pa.dictionary(pa.int64(), pa.large_string())),
+            pa.field("myfloat", pa.large_string()),
         ]
     )
 
@@ -1759,32 +1759,32 @@ def test_only_evolve_schema_when_enmr_is_extended(tmp_path):
     with soma.DataFrame.create(uri, schema=schema, domain=[[0, 4]]) as sdf:
         data = {}
         data["soma_joinid"] = [0, 1, 2, 3, 4]
-        data["foo"] = pd.Categorical(["a", "bb", "ccc", "bb", "a"])
-        data["bar"] = ["cat", "dog", "cat", "cat", "cat"]
+        data["myint"] = pd.Categorical(["a", "bb", "ccc", "bb", "a"])
+        data["myfloat"] = ["cat", "dog", "cat", "cat", "cat"]
         sdf.write(pa.Table.from_pydict(data))
 
     # +1 evolving the schema
     with soma.DataFrame.open(uri, "w") as sdf:
         data = {}
         data["soma_joinid"] = [0, 1, 2, 3, 4]
-        data["foo"] = pd.Categorical(["a", "bb", "ccc", "d", "a"])
-        data["bar"] = ["cat", "dog", "cat", "cat", "cat"]
+        data["myint"] = pd.Categorical(["a", "bb", "ccc", "d", "a"])
+        data["myfloat"] = ["cat", "dog", "cat", "cat", "cat"]
         sdf.write(pa.Table.from_pydict(data))
 
     # +0 no changes to enumeration values
     with soma.DataFrame.open(uri, "w") as sdf:
         data = {}
         data["soma_joinid"] = [0, 1, 2, 3, 4]
-        data["foo"] = pd.Categorical(["a", "bb", "ccc", "d", "a"])
-        data["bar"] = ["cat", "dog", "cat", "cat", "cat"]
+        data["myint"] = pd.Categorical(["a", "bb", "ccc", "d", "a"])
+        data["myfloat"] = ["cat", "dog", "cat", "cat", "cat"]
         sdf.write(pa.Table.from_pydict(data))
 
     # +0 no changes enumeration values
     with soma.DataFrame.open(uri, "w") as sdf:
         data = {}
         data["soma_joinid"] = [0, 1, 2, 3, 4]
-        data["foo"] = pd.Categorical(["a", "bb", "ccc", "d", "d"])
-        data["bar"] = ["cat", "dog", "cat", "cat", "cat"]
+        data["myint"] = pd.Categorical(["a", "bb", "ccc", "d", "d"])
+        data["myfloat"] = ["cat", "dog", "cat", "cat", "cat"]
         sdf.write(pa.Table.from_pydict(data))
 
     # total 3 fragment files
@@ -1899,3 +1899,47 @@ def test_bounds_on_somajoinid_domain(tmp_path):
     )
 
     assert soma.DataFrame.exists(uri)
+
+
+def test_pass_configs(tmp_path, arrow_schema):
+    uri = tmp_path.as_posix()
+
+    with soma.DataFrame.create(uri, schema=arrow_schema()) as sdf:
+        pydict = {}
+        pydict["soma_joinid"] = [0, 1, 2, 3, 4]
+        pydict["myint"] = [10, 20, 30, 40, 50]
+        pydict["myfloat"] = [4.1, 5.2, 6.3, 7.4, 8.5]
+        pydict["mystring"] = ["apple", "ball", "cat", "dog", "egg"]
+        pydict["mybool"] = [True, False, False, True, False]
+        rb = pa.Table.from_pydict(pydict)
+        sdf.tiledbsoma_resize_soma_joinid_shape(len(rb))
+        sdf.write(rb)
+
+    # Pass a custom config to open
+    with soma.DataFrame.open(
+        uri,
+        "r",
+        context=soma.SOMATileDBContext(
+            {"sm.mem.total_budget": "0", "sm.io_concurrency_level": "0"}
+        ),
+    ) as sdf:
+
+        # This errors out as 0 is not a valid value to set the total memory
+        # budget or nummber of threads
+        with pytest.raises(soma.SOMAError):
+            next(sdf.read())
+
+        # This still errors out because read still sees that the number of
+        # threads is 0 and therefore invalid
+        with pytest.raises(soma.SOMAError):
+            next(sdf.read(platform_config={"sm.mem.total_budget": "10000"}))
+
+        # With correct values, this reads without issue
+        next(
+            sdf.read(
+                platform_config={
+                    "sm.mem.total_budget": "10000",
+                    "sm.io_concurrency_level": "1",
+                }
+            )
+        )
diff --git a/apis/python/tests/test_dense_nd_array.py b/apis/python/tests/test_dense_nd_array.py
index dea231f295..693b96da4a 100644
--- a/apis/python/tests/test_dense_nd_array.py
+++ b/apis/python/tests/test_dense_nd_array.py
@@ -498,3 +498,45 @@ def test_read_to_unwritten_array(tmp_path, shape):
 
     with soma.DenseNDArray.open(uri, "r") as A:
         assert np.array_equal(np.ones(shape) * 255, A.read().to_numpy())
+
+
+def test_pass_configs(tmp_path):
+    uri = tmp_path.as_posix()
+
+    with soma.DenseNDArray.create(
+        tmp_path.as_posix(),
+        type=pa.uint8(),
+        shape=(2, 2),
+        context=SOMATileDBContext(timestamp=1),
+    ) as a:
+        a.write(
+            (slice(0, 2), slice(0, 2)),
+            pa.Tensor.from_numpy(np.zeros((2, 2), dtype=np.uint8)),
+        )
+
+    # Pass a custom config to open
+    with soma.DenseNDArray.open(
+        uri,
+        "r",
+        context=soma.SOMATileDBContext(
+            {"sm.mem.total_budget": "0", "sm.io_concurrency_level": "0"}
+        ),
+    ) as sdf:
+
+        # This errors out as 0 is not a valid value to set the total memory
+        # budget or nummber of threads
+        with pytest.raises(soma.SOMAError):
+            sdf.read()
+
+        # This still errors out because read still sees that the number of
+        # threads is 0 and therefore invalid
+        with pytest.raises(soma.SOMAError):
+            sdf.read(platform_config={"sm.mem.total_budget": "300000"})
+
+        # With correct values, this reads without issue
+        sdf.read(
+            platform_config={
+                "sm.mem.total_budget": "300000",
+                "sm.io_concurrency_level": "1",
+            }
+        )
diff --git a/apis/python/tests/test_export_multiscale_image.py b/apis/python/tests/test_export_multiscale_image.py
index bc6609164c..10b1f6a984 100644
--- a/apis/python/tests/test_export_multiscale_image.py
+++ b/apis/python/tests/test_export_multiscale_image.py
@@ -48,17 +48,31 @@ def sample_multiscale_image_2d(tmp_path_factory, sample_2d_data):
 
 
 @pytest.mark.parametrize(
-    "level,transform,expected_transformation",
+    "level,transform,expected_transformation,expected_transformation_key",
     [
+        (
+            0,
+            None,
+            sd.transformations.Identity(),
+            "image",
+        ),
+        (
+            2,
+            None,
+            sd.transformations.Scale([4, 4], ("x", "y")),
+            "image",
+        ),
         (
             0,
             somacore.IdentityTransform(("x_scene", "y_scene"), ("x_image", "y_image")),
             sd.transformations.Identity(),
+            "scene0",
         ),
         (
             2,
             somacore.IdentityTransform(("x_scene", "y_scene"), ("x_image", "y_image")),
             sd.transformations.Scale([4, 4], ("x", "y")),
+            "scene0",
         ),
         (
             0,
@@ -66,6 +80,7 @@ def sample_multiscale_image_2d(tmp_path_factory, sample_2d_data):
                 ("x_scene", "y_scene"), ("x_image", "y_image"), [0.25, 0.5]
             ),
             sd.transformations.Scale([4, 2], ("x", "y")),
+            "scene0",
         ),
         (
             2,
@@ -73,6 +88,7 @@ def sample_multiscale_image_2d(tmp_path_factory, sample_2d_data):
                 ("x_scene", "y_scene"), ("x_image", "y_image"), [0.25, 0.5]
             ),
             sd.transformations.Scale([16, 8], ("x", "y")),
+            "scene0",
         ),
         (
             0,
@@ -84,6 +100,7 @@ def sample_multiscale_image_2d(tmp_path_factory, sample_2d_data):
                 ("x", "y"),
                 ("x", "y"),
             ),
+            "scene0",
         ),
         (
             2,
@@ -100,6 +117,7 @@ def sample_multiscale_image_2d(tmp_path_factory, sample_2d_data):
                     ),
                 ]
             ),
+            "scene0",
         ),
     ],
 )
@@ -109,10 +127,12 @@ def test_export_image_level_to_spatial_data(
     level,
     transform,
     expected_transformation,
+    expected_transformation_key,
 ):
     image2d = soma_outgest.to_spatial_data_image(
         sample_multiscale_image_2d,
         level=level,
+        key="image",
         scene_id="scene0",
         scene_dim_map={"x_scene": "x", "y_scene": "y"},
         transform=transform,
@@ -131,12 +151,23 @@ def test_export_image_level_to_spatial_data(
     # Check the metadata.
     metadata = dict(image2d.attrs)
     assert len(metadata) == 1
-    assert metadata["transform"] == {"scene0": expected_transformation}
+    assert metadata["transform"] == {
+        expected_transformation_key: expected_transformation
+    }
 
 
 @pytest.mark.parametrize(
-    "transform,expected_transformation",
+    "transform,expected_transformation,expected_transformation_key",
     [
+        (
+            None,
+            [
+                sd.transformations.Identity(),
+                sd.transformations.Scale([2, 2], ("x", "y")),
+                sd.transformations.Scale([4, 4], ("x", "y")),
+            ],
+            "image",
+        ),
         (
             somacore.IdentityTransform(("x_scene", "y_scene"), ("x_image", "y_image")),
             [
@@ -144,6 +175,7 @@ def test_export_image_level_to_spatial_data(
                 sd.transformations.Scale([2, 2], ("x", "y")),
                 sd.transformations.Scale([4, 4], ("x", "y")),
             ],
+            "scene0",
         ),
         (
             somacore.ScaleTransform(
@@ -154,6 +186,7 @@ def test_export_image_level_to_spatial_data(
                 sd.transformations.Scale([8, 4], ("x", "y")),
                 sd.transformations.Scale([16, 8], ("x", "y")),
             ],
+            "scene0",
         ),
         (
             somacore.AffineTransform(
@@ -186,14 +219,20 @@ def test_export_image_level_to_spatial_data(
                     ]
                 ),
             ],
+            "scene0",
         ),
     ],
 )
 def test_export_full_image_to_spatial_data(
-    sample_multiscale_image_2d, sample_2d_data, transform, expected_transformation
+    sample_multiscale_image_2d,
+    sample_2d_data,
+    transform,
+    expected_transformation,
+    expected_transformation_key,
 ):
     image2d = soma_outgest.to_spatial_data_multiscale_image(
         sample_multiscale_image_2d,
+        key="image",
         scene_id="scene0",
         scene_dim_map={"x_scene": "x", "y_scene": "y"},
         transform=transform,
@@ -216,4 +255,6 @@ def test_export_full_image_to_spatial_data(
         # Check the metadata.
         metadata = dict(data_array.attrs)
         assert len(metadata) == 1
-        assert metadata["transform"] == {"scene0": expected_transformation[index]}
+        assert metadata["transform"] == {
+            expected_transformation_key: expected_transformation[index]
+        }
diff --git a/apis/python/tests/test_export_point_cloud_dataframe.py b/apis/python/tests/test_export_point_cloud_dataframe.py
index 58823ffee8..bd7f50f736 100644
--- a/apis/python/tests/test_export_point_cloud_dataframe.py
+++ b/apis/python/tests/test_export_point_cloud_dataframe.py
@@ -11,7 +11,7 @@
 
 gpd = pytest.importorskip("geopandas")
 soma_outgest = pytest.importorskip("tiledbsoma.io.spatial.outgest")
-spatialdata = pytest.importorskip("spatialdata")
+sd = pytest.importorskip("spatialdata")
 
 
 @pytest.fixture(scope="module")
@@ -35,21 +35,35 @@ def sample_point_cloud_dataframe_2d(tmp_path_factory):
     return point_cloud
 
 
-def test_export_to_shapes_2d(sample_point_cloud_dataframe_2d):
+@pytest.mark.parametrize(
+    "transform,expected_transformation",
+    [
+        (None, {"point_cloud": sd.transformations.Identity()}),
+        (
+            somacore.IdentityTransform(
+                ("x_scene", "y_scene"), ("x_points", "y_points")
+            ),
+            {"scene0": sd.transformations.Identity()},
+        ),
+    ],
+)
+def test_export_to_shapes_2d(
+    sample_point_cloud_dataframe_2d, transform, expected_transformation
+):
     """Test exporting a simple point cloud to a SpatialData shape model."""
     # Export PointCloudDataFrame to shapes.
     shape = soma_outgest.to_spatial_data_shapes(
         sample_point_cloud_dataframe_2d,
         scene_id="scene0",
+        key="point_cloud",
         scene_dim_map={"x_scene": "x", "y_scene": "y"},
         soma_joinid_name="obs_id",
-        transform=somacore.IdentityTransform(
-            ("x_scene", "y_scene"), ("x_points", "y_points")
-        ),
+        transform=transform,
     )
 
     # Check this is valid storage for the SpatialData "Shapes" model.
-    spatialdata.models.ShapesModel.validate(shape)
+    model_schema = sd.models.get_model(shape)
+    assert model_schema == sd.models.ShapesModel
 
     # Check the dataframe.
     expected = gpd.GeoDataFrame.from_dict(
@@ -65,27 +79,39 @@ def test_export_to_shapes_2d(sample_point_cloud_dataframe_2d):
 
     # Check the metadata.
     metadata = dict(shape.attrs)
-    for key, val in metadata.items():
-        print(f"{key}: {val}")
     assert len(metadata) == 1
-    assert metadata["transform"] == {"scene0": spatialdata.transformations.Identity()}
+    assert metadata["transform"] == expected_transformation
 
 
-def test_export_to_points_2d(sample_point_cloud_dataframe_2d):
+@pytest.mark.parametrize(
+    "transform,expected_transformation",
+    [
+        (None, {"point_cloud": sd.transformations.Identity()}),
+        (
+            somacore.IdentityTransform(
+                ("x_scene", "y_scene"), ("x_points", "y_points")
+            ),
+            {"scene0": sd.transformations.Identity()},
+        ),
+    ],
+)
+def test_export_to_points_2d(
+    sample_point_cloud_dataframe_2d, transform, expected_transformation
+):
     """Test exporting a simple point cloud to a SpatialData shape model."""
     # Export PointCloudDataFrame to shapes.
     points = soma_outgest.to_spatial_data_points(
         sample_point_cloud_dataframe_2d,
+        key="point_cloud",
         scene_id="scene0",
         scene_dim_map={"x_scene": "x", "y_scene": "y"},
         soma_joinid_name="obs_id",
-        transform=somacore.IdentityTransform(
-            ("x_scene", "y_scene"), ("x_points", "y_points")
-        ),
+        transform=transform,
     )
 
     # Check this is valid storage for the SpatialData "Points" model.
-    spatialdata.models.PointsModel.validate(points)
+    model_schema = sd.models.get_model(points)
+    assert model_schema == sd.models.PointsModel
 
     # Check the dataframe.
     expected = pd.DataFrame.from_dict(
@@ -99,7 +125,5 @@ def test_export_to_points_2d(sample_point_cloud_dataframe_2d):
 
     # Check the metadata.
     metadata = dict(points.attrs)
-    for key, val in metadata.items():
-        print(f"{key}: {val}")
     assert len(metadata) == 1
-    assert metadata["transform"] == {"scene0": spatialdata.transformations.Identity()}
+    assert metadata["transform"] == expected_transformation
diff --git a/apis/python/tests/test_scene.py b/apis/python/tests/test_scene.py
index 625ac7d984..f8e6ad4b1c 100644
--- a/apis/python/tests/test_scene.py
+++ b/apis/python/tests/test_scene.py
@@ -111,7 +111,39 @@ def test_measurement_with_var_scene(tmp_path):
     assert soma.DataFrame.exists(obs_scene_uri)
 
 
-def test_scene_coord_space(tmp_path):
+def test_scene_coord_space_at_create(tmp_path):
+    uri = tmp_path.as_uri()
+
+    coord_space = soma.CoordinateSpace(
+        [
+            soma.Axis(name="x"),
+            soma.Axis(name="y"),
+        ]
+    )
+    coord_space_json = """
+    [
+        {"name": "x", "unit": null},
+        {"name": "y", "unit": null}
+    ]
+    """
+
+    with soma.Scene.create(uri, coordinate_space=("x", "y")) as scene:
+
+        # Reserved metadata key should not be settable?
+        # with pytest.raises(soma.SOMAError):
+        #     scene.metadata["soma_coordinate_space"] = coord_space_json
+
+        scene.coordinate_space = coord_space
+        assert scene.coordinate_space == coord_space
+        assert json.loads(scene.metadata["soma_coordinate_space"]) == json.loads(
+            coord_space_json
+        )
+
+    with soma.Scene.open(uri) as scene:
+        assert scene.coordinate_space == coord_space
+
+
+def test_scene_coord_space_after_create(tmp_path):
     uri = tmp_path.as_uri()
 
     coord_space = soma.CoordinateSpace(
diff --git a/apis/python/tests/test_sparse_nd_array.py b/apis/python/tests/test_sparse_nd_array.py
index 403e806b6d..cedd88c68e 100644
--- a/apis/python/tests/test_sparse_nd_array.py
+++ b/apis/python/tests/test_sparse_nd_array.py
@@ -1883,3 +1883,47 @@ def test_global_writes(tmp_path):
                 data,
                 platform_config=soma.TileDBCreateOptions(),
             )
+
+
+def test_pass_configs(tmp_path):
+    uri = tmp_path.as_posix()
+
+    with soma.SparseNDArray.create(
+        tmp_path.as_posix(), type=pa.uint8(), shape=(3,)
+    ) as a:
+        data = pa.Table.from_pydict(
+            {
+                "soma_dim_0": pa.array([0, 1, 2], type=pa.int64()),
+                "soma_data": pa.array([1, 2, 3], type=pa.uint8()),
+            }
+        )
+        a.write(data)
+
+    # Pass a custom config to open
+    with soma.SparseNDArray.open(
+        uri,
+        "r",
+        context=soma.SOMATileDBContext(
+            {"sm.mem.total_budget": "0", "sm.io_concurrency_level": "0"}
+        ),
+    ) as sdf:
+
+        # This errors out as 0 is not a valid value to set the total memory
+        # budget or nummber of threads
+        with pytest.raises(soma.SOMAError):
+            next(sdf.read().tables())
+
+        # This still errors out because read still sees that the number of
+        # threads is 0 and therefore invalid
+        with pytest.raises(soma.SOMAError):
+            next(sdf.read(platform_config={"sm.mem.total_budget": "10000"}).tables())
+
+        # With correct values, this reads without issue
+        next(
+            sdf.read(
+                platform_config={
+                    "sm.mem.total_budget": "10000",
+                    "sm.io_concurrency_level": "1",
+                }
+            ).tables()
+        )
diff --git a/apis/python/tests/test_spatial_outgest_util.py b/apis/python/tests/test_spatial_outgest_util.py
index 0c771ee998..ca224f26ce 100644
--- a/apis/python/tests/test_spatial_outgest_util.py
+++ b/apis/python/tests/test_spatial_outgest_util.py
@@ -2,7 +2,7 @@
 import pytest
 import somacore
 
-soma_outgest = pytest.importorskip("tiledbsoma.experimental.outgest")
+soma_outgest = pytest.importorskip("tiledbsoma.io.spatial.outgest")
 sd = pytest.importorskip("spatialdata")
 
 
diff --git a/apis/r/tests/testthat/test-00-Blockwise.R b/apis/r/tests/testthat/test-00-Blockwise.R
index 66e0eb511d..a34c20759d 100644
--- a/apis/r/tests/testthat/test-00-Blockwise.R
+++ b/apis/r/tests/testthat/test-00-Blockwise.R
@@ -152,6 +152,7 @@ test_that("Blockwise iterator for sparse matrices", {
   expect_error(bi$sparse_matrix("R"))
 
   expect_s3_class(it <- bi$sparse_matrix(), "BlockwiseSparseReadIter")
+  expect_false(it$reindexable)
   expect_false(it$read_complete())
 
   for (i in seq.int(1L, ceiling(it$coords_axis$length() / it$coords_axis$stride))) {
@@ -200,6 +201,7 @@ test_that("Sparse matrix blockwise iterator: re-indexed", {
   query <- exp$axis_query("RNA")
   xrqry <- query$X("data")
 
+  # Reindex only on major axis
   expect_s3_class(
     bi <- xrqry$blockwise(axis = ax, size = sz, reindex_disable_on_axis = NA),
     "SOMASparseNDArrayBlockwiseRead"
@@ -216,7 +218,11 @@ test_that("Sparse matrix blockwise iterator: re-indexed", {
   for (i in seq.int(1L, ceiling(it$coords_axis$length() / it$coords_axis$stride))) {
     mat <- it$read_next()
     expect_s4_class(mat, "TsparseMatrix")
-    expect_identical(dim(mat), rev(dim(obj)))
+    dims <- c(
+      ifelse(it$read_complete(), yes = ncol(obj) %% sz, no = sz),
+      nrow(obj)
+    )
+    expect_identical(dim(mat), dims)
     expect_true(min(mat@i) >= 0L)
     expect_true(max(mat@i) <= sz)
     strider <- attr(mat, "coords")$soma_dim_0
@@ -225,6 +231,7 @@ test_that("Sparse matrix blockwise iterator: re-indexed", {
     expect_true(strider$end < sz * i)
   }
 
+  # Reindex on all axes
   expect_s3_class(
     bi <- suppressWarnings(xrqry$blockwise(
       axis = ax,
@@ -242,7 +249,11 @@ test_that("Sparse matrix blockwise iterator: re-indexed", {
   for (i in seq.int(1L, ceiling(it$coords_axis$length() / it$coords_axis$stride))) {
     mat <- it$read_next()
     expect_s4_class(mat, "TsparseMatrix")
-    expect_identical(dim(mat), rev(dim(obj)))
+    dims <- c(
+      ifelse(it$read_complete(), yes = ncol(obj) %% sz, no = sz),
+      nrow(obj)
+    )
+    expect_identical(dim(mat), dims)
     expect_true(min(mat@i) >= 0L)
     expect_true(max(mat@i) <= sz)
   }
@@ -270,14 +281,19 @@ test_that("Blockwise iterate through full array", {
   obs_stride <- n_obs %/% n_chunks
   it <- exp$ms$get("RNA")$X$get(X_layer)$read()$blockwise(axis = 0L, size = obs_stride)$sparse_matrix()
   expect_false(it$read_complete())
+  i <- 1L
   while (!it$read_complete()) {
+    i <- i + 1L
     mat <- it$read_next()
     expect_s4_class(mat, "dgTMatrix")
     expect_identical(ncol(mat), n_var)
-    expect_identical(nrow(mat), n_obs)
+    nobs <- ifelse(it$read_complete(), yes = n_obs %% obs_stride, no = obs_stride)
+    expect_identical(nrow(mat), nobs, expected.label = nobs)
+    ncoords <- ifelse(it$read_complete(), yes = n_obs %% n_chunks, no = obs_stride)
     expect_identical(
       length(attr(mat, "coords")$soma_dim_0),
-      ifelse(it$read_complete(), yes = n_obs %% n_chunks, no = obs_stride)
+      ncoords,
+      expected.label = ncoords
     )
   }
   expect_true(it$read_complete())
@@ -289,11 +305,14 @@ test_that("Blockwise iterate through full array", {
   while (!it$read_complete()) {
     mat <- it$read_next()
     expect_s4_class(mat, "dgTMatrix")
-    expect_identical(ncol(mat), n_var)
+    nvar <- ifelse(it$read_complete(), yes = n_var %% var_stride, no = var_stride)
+    expect_identical(ncol(mat), nvar, expected.label = nvar)
     expect_identical(nrow(mat), n_obs)
+    ncoords <- ifelse(it$read_complete(), yes = n_var %% n_chunks, no = var_stride)
     expect_identical(
       length(attr(mat, "coords")$soma_dim_1),
-      ifelse(it$read_complete(), yes = n_var %% n_chunks, no = var_stride)
+      ncoords,
+      expected.label = ncoords
     )
   }
   expect_true(it$read_complete())
diff --git a/apis/r/tests/testthat/test-11-shape.R b/apis/r/tests/testthat/test-11-shape.R
index 0b510ea071..088e9fa73e 100644
--- a/apis/r/tests/testthat/test-11-shape.R
+++ b/apis/r/tests/testthat/test-11-shape.R
@@ -360,6 +360,21 @@ test_that("SOMADataFrame domain mods", {
         myint = c(20, 50),
         myfloat = c(0.0, 6.0)
     )
+
+    # -- first check dry run
+    expect_no_condition(sdf$change_domain(domain_for_create, check_only = TRUE))
+    sdf$close()
+
+    check <- list(
+      soma_joinid = c(0, 3),
+      mystring = c("", ""), # this is how it reads back
+      myint = c(20, 50),
+      myfloat = c(0.0, 6.0)
+    )
+    expect_equal(sdf$domain(), check)
+    sdf$close()
+
+    sdf <- SOMADataFrameOpen(uri, "WRITE")
     expect_no_condition(sdf$change_domain(new_domain))
 
     # Shrink
@@ -409,10 +424,11 @@ test_that("SOMADataFrame domain mods", {
     sdf$close()
 
     # Check for success
-    sdf <- SOMADataFrameOpen(uri, "WRITE")
+    sdf <- SOMADataFrameOpen(uri, "READ")
     dom <- sdf$domain()
     expect_equal(sdf$domain(), new_domain)
     sdf$close()
+
 })
 
 test_that("SOMASparseNDArray shape", {
@@ -555,7 +571,12 @@ test_that("SOMADenseNDArray shape", {
     expect_error(ndarray$write(mat))
     ndarray$close()
 
-    # Test resize up
+    # Test resize up, dry run
+    new_shape <- c(501, 601)
+    expect_no_error(reason_string <- ndarray$resize(new_shape, check_only = TRUE))
+    expect_equal(reason_string, "")
+
+    # Test resize up, for real
     new_shape <- c(500, 600)
     expect_no_error(ndarray$resize(new_shape))
 
diff --git a/libtiledbsoma/test/unit_soma_group.cc b/libtiledbsoma/test/unit_soma_group.cc
index 36775a9595..2b7c5b1189 100644
--- a/libtiledbsoma/test/unit_soma_group.cc
+++ b/libtiledbsoma/test/unit_soma_group.cc
@@ -100,7 +100,7 @@ std::tuple<std::string, uint64_t> create_array(
     // Open array for writing
     Array array(ctx, uri, TILEDB_WRITE, TemporalPolicy(TimeTravel, timestamp));
     if (LOG_DEBUG_ENABLED()) {
-        array.schema().dump();
+        std::cout << array.schema();
     }
 
     // Generate fragments in random order