n = 500
df = md.DataFrame(
    {
        "coord_x": mt.random.uniform(-3, 3, size=n),
        "coord_y": mt.random.uniform(30, 50, size=n),
        "observations": mt.random.randint(1, 5, size=n),
    }
)
ax = df.plot.hexbin(
    x="coord_x",
    y="coord_y",
    C="observations",
    reduce_C_function=mt.sum,
    gridsize=10,
    cmap="viridis",
)