Chunking Warning in xarray: Understanding and Optimizing Performance

[rjavierch]

Hello xarray community,

I'm currently working on a project involving mid-large datasets and have encountered a warning related to chunking that I'm hoping to get some clarification on. The warning message I'm receiving is as follows:

UserWarning: The specified chunks separate the stored chunks along dimension "x" starting at index 5000. This could degrade performance. Instead, consider rechunking after loading.

I do not understand the reason for this warning, especially since when I initially open the dataset with no specified chunk sizes and execute the .chunksizes method, the result indicates that there are no chunks. However, when I attempt to specify chunks along the "x" dimension starting at index 1000, I receive this warning.

Could someone please provide some insights into why this warning occurs and how I can address it? Specifically, I'm interested in understanding:

The implications of the warning message in terms of performance degradation.
Best practices for determining an optimal chunking strategy, especially when starting chunking at a specific index.

Thank you in advance for your help and insights!

The dataset I am currently using is the raster file of the Harmonized World Soil Database v2.0

This is the code that shows the warning:

chunk_size = {
    "x": 5000,
    "y": 5000,
}
# Does not requires chunks as was written already chunked
hwsd2_rast = xr.open_dataset(
    "HWSD2/HWSD2.bil",
    chunks=chunk_size,
    engine="rasterio",
    )
hwsd2_rast

With hwsd2_rast.chunksizes output:

UserWarning:

The specified chunks separate the stored chunks along dimension "x" starting at index 5000. This could degrade performance. Instead, consider rechunking after loading.

Frozen({'band': (1,), 'y': (5000, 5000, 5000, 5000, 1600), 'x': (5000, 5000, 5000, 5000, 5000, 5000, 5000, 5000, 3200)})

If I remove the chunks for "x" only, there is no warning. However, these are the chunk sizes for hwsd2_rast.chunksizes

Frozen({'band': (1,), 'y': (5000, 5000, 5000, 5000, 1600), 'x': (43200,)})

Thank you!

[keewis]

can you check the encoding for a variable in your dataset that has the dimensions in question? Something like

hwsd2_rast["<variable1>"].encoding

This will tell you the size of the chunks on disk, which are the smallest unit of what we can load into memory.

If your dataset is not chunked along x on disk, choosing a chunksize smaller than that will load the full chunk (all 43200 values along that axis) and drop the ones that are not selected.

[rjavierch]

Sure, I have excecuted:

hwsd2_rast["<variable1>"].encoding

for coordinates "x" and "y" and variable (unique) "band_data" with the respective following outputs:

{ }

{ }

{'dtype': 'uint16',
 'scale_factor': 1.0,
 'add_offset': 0.0,
 '_FillValue': 65535.0,
 'grid_mapping': 'spatial_ref',
 'preferred_chunks': {'y': 1, 'x': 43200, 'band': 1},
 'source': 'HWSD2/HWSD2.bil',
 'rasterio_dtype': 'uint16'}

I think, based on what you have explained and what I understood, the data was stored in one unique chunk for coordinate "y" and in 43200 chunks for coordinate "x", right?

Aside, according to the warning, I should then implement my code by ignoring the chunks argument in the open_dataset function and rechunking later like in the following code?

chunk_size = {
    "x": 5000,
    "y": 5000,
}
# Does not requires chunks as was written already chunked
hwsd2_rast = xr.open_dataset(
    "HWSD2/HWSD2.bil",
    engine="rasterio",
    )
hwsd2_rast = hwsd2_rast.chunk(chunk_size)
hwsd2_rast

[keewis]

not quite: the chunk size tells you the size of the chunks, not the number of chunks along that dimension. In this case, your dataset appears to store the data row by row. Thus, instead of somewhat rectangular chunks, it appears you have a single chunk along x and many 1-sized chunks along y.

With the code you posted, you'll load 5000 of these on-disk chunks into memory on a worker, then split them into 9 smaller chunks along x (which is what the warning recommended). For that to work properly, however, you need to specify chunks={} in the open_dataset call, which will give you the chunksize on disk.

[rjavierch]

I see! Now it is clearer for me.
Thanks a lot for your time!