[Numpy-discussion] weighted mean; weighted standard error of the mean (sem)
Thu Sep 9 21:38:18 CDT 2010
Excerpts from cpblpublic's message of Thu Sep 09 22:22:05 -0400 2010:
> I am looking for some reaally basic statistical tools. I have some
> sample data, some sample weights for those measurements, and I want to
> calculate a mean and a standard error of the mean.
> Here are obvious places to look:
> It seems to me that numpy's "mean" and "average" functions have their
> names backwards. That is, often a mean is defined more generally than
> average, and includes the possibility of weighting, but in this case
> it is "average" that has a weights argument. Can these functions be
> merged/renamed/deprecated in the future? It's clear to me that "mean"
> should allow for weights.
> None of these modules, above, offer standard error of the mean which
> incorporates weights. scipy.stats.sem() doesn't, and that's the closest
> thing. numpy's "var" doesn't allow weights.
> There aren't any weighted variances in the above modules.
> Again, are there favoured codes for these functions? Can they be
> incorporated appropriately in the future?
> Most immediately, I'd love to get code for weighted sem. I'll write it
> otherwise, but it might be crude and dumb...
> Chris Barrington-Leigh
This code below should do what you want. It is part of the stat
sub-package of esutil http://code.google.com/p/esutil/
Hope this helps,
Erin Scott Sheldon
Brookhaven National Laboratory
def wmom(arrin, weights_in, inputmean=None, calcerr=False, sdev=False):
Calculate the weighted mean, error, and optionally standard deviation of
an input array. By default error is calculated assuming the weights are
1/err^2, but if you send calcerr=True this assumption is dropped and the
error is determined from the weighted scatter.
wmean,werr = wmom(arr, weights, inputmean=None, calcerr=False, sdev=False)
arr: A numpy array or a sequence that can be converted.
weights: A set of weights for each elements in array.
An input mean value, around which them mean is calculated.
Calculate the weighted error. By default the error is calculated as
1/sqrt( weights.sum() ). If calcerr=True it is calculated as sqrt(
(w**2 * (arr-mean)**2).sum() )/weights.sum()
If True, also return the weighted standard deviation as a third
element in the tuple.
wmean, werr: A tuple of the weighted mean and error. If sdev=True the
tuple will also contain sdev: wmean,werr,wsdev
Converted from IDL: 2006-10-23. Erin Sheldon, NYU
# no copy made if they are already arrays
arr = numpy.array(arrin, ndmin=1, copy=False)
# Weights is forced to be type double. All resulting calculations
# will also be double
weights = numpy.array(weights_in, ndmin=1, dtype='f8', copy=False)
wtot = weights.sum()
# user has input a mean value
if inputmean is None:
wmean = ( weights*arr ).sum()/wtot
# how should error be calculated?
werr2 = ( weights**2 * (arr-wmean)**2 ).sum()
werr = numpy.sqrt( werr2 )/wtot
werr = 1.0/numpy.sqrt(wtot)
# should output include the weighted standard deviation?
wvar = ( weights*(arr-wmean)**2 ).sum()/wtot
wsdev = numpy.sqrt(wvar)
More information about the NumPy-Discussion