def color_stats(lab_image, mask=None):
    """Get color stats of a L*a*b* image (mean, std-dev)"""
    if mask is None:
        # Image is always in 1 dimension to make it simpler
        lab_image = lab_image.reshape(-1, lab_image.shape[-1])
    else:
        lab_image = lab_image[mask.nonzero()]
    return lab_image.mean(axis=0), lab_image.std(axis=0)