Wnioskowanie Statystyczne

Jak wyznaczyć przedział ufności dla dowolnej statystyki?

Przy użyciu symulacji!

Mając dany zestaw \(N\) wartości \(x_1 \ldots x_N\) z dowolnego rozkładu, możemy zastosować metodę bootstrap w celu znalezienia przedziału ufności dla dowolnie wybranej statystyki (funkcji próby), czyli np. mediany. Jako miarę pewności oszacowania, ustalamy poziom istotności \(\alpha\) jako maksymalne dopuszczane przez nas prawdopodobieństwo, że nasz przedział nie obejmie prawdziwej wartości estymowanego parametru. Zazwyczaj przyjmujemy w tym celu wartość \(\alpha = 0.05\) lub mniejszą.

W tym celu wielokrotnie powtarzamy następujące kroki:

import numpy
import scipy.stats

# poziom istotności
alpha = 0.01

# dane wejściowe
N = 5
X = numpy.array([1.2, 2.3, 3.1, 1.7, 2.5])
ile_powtorzen = 10000

# --- WERSJA 1 ---
statystyki = numpy.zeros(ile_powtorzen)
for k in xrange(ile_powtorzen) :
  Y = numpy.zeros(N)
  for j in xrange(N) :
    Y[j] = X[numpy.random.randint(N)]
  statystyki[k] = TUTAJ OBLICZAMY WARTOŚĆ STATYSTYKI OD Y

# --- WERSJA 2 ---
statystyki = numpy.zeros(ile_powtorzen)
for k in xrange(ile_powtorzen) :
  Y = X[numpy.random.randint(N, size=N)]
  statystyki[k] = TUTAJ OBLICZAMY WARTOŚĆ STATYSTYKI OD Y

Wynikiem powyższego etapu powinna być statystyk z poszczególnych „syntetycznych” prób (zestawów danych). Z tej listy możemy obliczyć kwantyle rzędu \(\frac{\alpha}{2}\) i \((1-\frac{\alpha}{2})\); są one już krańcami szukanego przedziału ufności.

# znajdujemy przedział ufności,
# pamiętajmy że funkcja scoreatpercentile przyjmuje wartości w procentach
xL = scipy.stats.scoreatpercentile(wariancje, 100.0 * alpha/2)
xP = scipy.stats.scoreatpercentile(wariancje, 100.0 * (1-alpha/2))
print [xL, xP]
autor: Piotr Różański, ostatnia modyfikacja: 10.04.2016