Why does statistics.variance use 'unbiased' sample variance by default?
<p>I have recently started using the statistics module for python.</p>

<p>I've noticed that by default the variance() method returns the 'unbiased' variance or sample variance:</p>

<pre><code>import statistics as st
from random import randint

def myVariance(data):
    # finds the variance of a given set of numbers
    xbar = st.mean(data)
    return sum([(x - xbar)**2 for x in data])/len(data)

def myUnbiasedVariance(data):
    # finds the 'unbiased' variance of a given set of numbers (divides by N-1) 
    xbar = st.mean(data)
    return sum([(x - xbar)**2 for x in data])/(len(data)-1)

population = [randint(0, 1000) for i in range(0,100)]

print myVariance(population)

print myUnbiasedVariance(population)

print st.variance(population)
</code></pre>

<p>output:</p>

<pre><code>81295.8011
82116.9708081
82116.9708081
</code></pre>

<p>This seems odd to me. I guess a lot of the time people are working with samples so they want a sample variance, but i would expect the default function to calculate a population variance. Does anyone know why this is?</p>

<p>I would argue that almost all the time when people estimate the variance from data they work with a sample.  And, by the definition of <em>unbiased estimate</em>, the expected value of the unbiased estimate of the variance equals the population variance.</p>

<p>In your code, you use <code>random.randint(0, 1000)</code>, which samples from a discrete uniform distribution with 1001 possible values and variance 1000*1002/12 = 83500 (see, e.g., <a href="http://mathworld.wolfram.com/DiscreteUniformDistribution.html" rel="nofollow">MathWorld</a>).  Here code that shows that, on average and when using samples as input, <code>statistics.variance()</code> gets closer to the population variance than <code>statistics.pvariance()</code>:</p>

<pre><code>import statistics as st, random, numpy as np

var, pvar = [], []
for i in range(10000):
  smpl = [random.randint(0, 1000) for j in range(10)]
  var.append(st.variance(smpl))
  pvar.append(st.pvariance(smpl))

print "mean variance(sample):  %.1f" %np.mean(var)
print "mean pvariance(sample): %.1f" %np.mean(pvar)
print "pvariance(population):  %.1f" %st.pvariance(range(1001))
</code></pre>

<p>Here sample output:</p>

<pre><code>mean variance(sample):  83626.0
mean pvariance(sample): 75263.4
pvariance(population):  83500.0
</code></pre>

python statistics variance