ols - ¿Cómo realizar una prueba de bondad de ajuste de ji cuadrado utilizando bibliotecas científicas en Python?

stats scipy python

Question

Supongamos que tengo algunos datos que obtuve empíricamente:

from scipy import stats
size = 10000
x = 10 * stats.expon.rvs(size=size) + 0.2 * np.random.uniform(size=size)

Se distribuye de forma exponencial (con algo de ruido) y quiero verificar esto mediante una prueba de bondad de ajuste (GoF) de ji cuadrado. ¿Cuál es la forma más sencilla de hacerlo utilizando las bibliotecas científicas estándar en Python (por ejemplo, scipy o statsmodels) con la menor cantidad de pasos y suposiciones manuales?

Puedo ajustar un modelo con:

param = stats.expon.fit(x)
plt.hist(x, normed=True, color=''white'', hatch=''/'')
plt.plot(grid, distr.pdf(np.linspace(0, 100, 10000), *param))

Es muy elegante calcular la prueba de Kolmogorov-Smirnov .

>>> stats.kstest(x, lambda x : stats.expon.cdf(x, *param))
(0.0061000000000000004, 0.85077099515985011)

Sin embargo, no puedo encontrar una buena manera de calcular la prueba de ji cuadrado.

Hay una función GoF de ji cuadrado en statsmodel , pero asume una distribución discreta (y la distribución exponencial es continua).

El tutorial oficial de scipy.stats solo cubre un caso para una distribución personalizada y las probabilidades se crean mediante el uso de muchas expresiones (npoints, npointsh, nbound, normbound), por lo que no tengo claro cómo hacerlo para otras distribuciones. Los ejemplos de chisquare suponen que los valores esperados y DoF ya se han obtenido.

Además, no estoy buscando una manera de realizar la prueba "manualmente" como ya se explicó aquí , pero me gustaría saber cómo aplicar una de las funciones de la biblioteca disponibles.

Answer 1

¿Por qué necesitas "verificar" que es exponencial? ¿Estás seguro de que necesitas una prueba estadística? Puedo garantizar bastante que no es en última instancia exponencial y la prueba sería significativa si tuviera suficientes datos, haciendo que la lógica de usar la prueba sea más bien forzada. Le puede ayudar a leer este hilo de CV: ¿Es la prueba de normalidad ''esencialmente inútil''? , o mi respuesta aquí: Pruebas de heterocedasticidad con muchas observaciones .

Por lo general, es mejor usar un qq-plot y / o pp-plot (dependiendo de si está preocupado por el ajuste en las colas o en la mitad de la distribución, vea mi respuesta aquí: PP-plots vs. QQ-plots ). La información sobre cómo hacer qq-plot en Python SciPy se puede encontrar en este subproceso de SO: Quantile-Quantile plot usando SciPy

Answer 2

Una solución aproximada para contenedores de probabilidad igual:

• Estimar los parámetros de la distribución.
• Use el cdf inverso, ppf si es un scipy.stats.distribution, para obtener los límites para una cuadrícula de probabilidad regular, por ejemplo, distribution.ppf(np.linspace(0, 1, n_bins + 1), *args)
• Luego, use np.histogram para contar el número de observaciones en cada casilla

A continuación, utilice la prueba de chisquare en las frecuencias.

Una alternativa sería encontrar los bordes de los intervalos de los percentiles de los datos ordenados, y usar el cdf para encontrar las probabilidades reales.

Esto es solo aproximado, ya que la teoría para la prueba de chisquare supone que los parámetros se estiman por máxima verosimilitud en los datos agrupados. Y no estoy seguro de si la selección de límites basados ​​en los datos afecta la distribución asintótica.

No he investigado esto en mucho tiempo. Si una solución aproximada no es lo suficientemente buena, entonces recomendaría que haga la pregunta en stats.stackexchange.