Downscaling arrays with categorical data (aka majority or mode resampling)

I have two global geospatial datasets with categorical data: one has 16 categories of land cover types and the other represents administrative regions of the world with nearly 400,000 different subregions. Both are global grids with resolution 0.01 degrees, stored in integer arrays of size 18000x36000.

I’m looking for a performant way of resampling these arrays to a smaller size, e.g. 1800x3600. It should use majority/mode resampling, i.e. choose the most common category within each subarea of the source array. (Algorithms that do antialiasing or interpolation, e.g. restrict or imresize from Images.jl, can’t be used since they would corrupt the categorical data.)

I know I can do this using GDAL (specifically gdal_translate) which has Julia wrappers in GDAL.jl and ArchGDAL.jl, but is there a package that does this kind of resampling in pure Julia?