Incorrect array statistics with (coverage) weights

[j08lue]

Incorrect array statistics with (coverage) weights

[vincentsarago]

Thanks @j08lue 🙏

count - what is that for weighted data? Count of all finite values where weights are non-zero? Or the sum of the weights?

I would sum the weights I think, I think it's the current behaviour: https://github.com/cogeotiff/rio-tiler/blob/main/rio_tiler/utils.py#L139

min, max, minority, majority, unique and also histogram - we should probably only take values into account where weights are non-zero.

so in

rio-tiler/rio_tiler/utils.py

Lines 134 to 163 in 73175d0

	array = data[b] * coverage
	else:
	array = data[b]
	count = array.count() if coverage is None else coverage.sum()
	if valid_pixels:
	majority = float(keys[counts.tolist().index(counts.max())].tolist())
	minority = float(keys[counts.tolist().index(counts.min())].tolist())
	else:
	majority = numpy.nan
	minority = numpy.nan
	output.append(
	{
	"min": float(data[b].min()),
	"max": float(data[b].max()),
	"mean": float(array.mean()),
	"count": float(count),
	"sum": float(array.sum()),
	"std": float(array.std()),
	"median": float(numpy.ma.median(array)),
	"majority": majority,
	"minority": minority,
	"unique": float(counts.size),
	**dict(zip(percentiles_names, percentiles_values)),
	"histogram": histogram,
	**info_px,
	}
	)

we should use array instead of data

[j08lue]

I would sum the weights

Yes, that must be the right way.

we should use array instead of data

Hmm, no, I think we need to get rid of multiplying coverage with data for all stats and only do that for mean.

What I meant was to make sure we expand the mask of data by the points where coverage equals zero. Since we compress that array (take out masked elements), we would then also exclude pixels with no coverage.

rio-tiler/rio_tiler/utils.py

Line 97 in 73175d0

keys, counts = numpy.unique(data[b].compressed(), return_counts=True)

Ideally, data.mask and coverage == 0 should be identical anyways, if data.mask is the rasterized query polygon and coverage is a refined version of the mask around the edges. But we need to be extra sure. Like

data.mask[:, coverage == 0] = True