Matrix Index is tilted using combine_by_coords

My calculations return a strange tilted index. Why does this happen?

What happened:
I combined several user defined chunked netcdf data (900 chunks) into one dataset. For this I used the default combine_by_coords in ds.open_mf_dataset(). My result was a tilted grid index - upper left corner i=0, j=1167.

Beforehand I calculated some indices on these chunks and combined them with the default combine='by_coords' in ds.open_mf_dataset() but also tested the combine='nested' separately.

The ones where I used default combine='by_coords' for all functions returned the tilted index.
The ones where I used combine='nested' beforehand and then default combine='by_coords' returned the correct index.

What you expected to happen:
No tilted index.

Minimal Complete Verifiable Example:

##returning wrong index
#calculating some climatic indices on numbered chunks and combined them to one ds per chunk
with xr.open_mfdataset(pathtofile+'annual*'+chunknumber+'.nc', chunks=-1, parallel=True, engine='h5netcdf') as ds:
    ds.to_netcdf(pathtofile, format="NETCDF4_CLASSIC", engine="netcdf4")

##combining all
with xr.open_mfdataset(pathtofile+'climateAnnual*.nc', chunks=-1, parallel=True, engine='h5netcd') as ds: 
    ds.to_netcdf(pathtofile, format="NETCDF4_CLASSIC", engine="netcdf4")

##################################################

##returning correct index
#calculating some climatic indices on numbered chunks and combined them to one ds per chunk
with xr.open_mfdataset(pathtofile+'annual*'+chunknumber+'.nc', chunks=-1, parallel=True, engine='h5netcdf', combine='nested') as ds:
    ds.to_netcdf(pathtofile, format="NETCDF4_CLASSIC", engine="netcdf4")

##combining all
with xr.open_mfdataset(pathtofile+'climateAnnual*.nc', chunks=-1, parallel=True, engine='h5netcd') as ds: 
    ds.to_netcdf(pathtofile, format="NETCDF4_CLASSIC", engine="netcdf4")

Anything else we need to know?:

Environment:
Python 3.7.4

Start with xarray.open_mfdataset, especially the default combine_by_coords path, and reproduce the result using the provided chunked NetCDF workflow. Compare it with combine='nested' and inspect how the coordinates and indices are ordered across the 900 chunks. Done means the cause of the tilted grid index is identified and the behavior is corrected or clearly documented.