¿Cómo interpreto la matriz de covarianza de un ajuste de curva?

No soy demasiado bueno en estadísticas, así que disculpa si esta es una pregunta simplista. Estoy ajustando una curva a algunos datos, y en ocasiones los datos de mi mejor encaja una exponencial negativa en la forma , y, a veces el ajuste es más cercano a un * e ( - b * x 2 ) + c . Sin embargo, a veces ambos fallan, y me gustaría volver a un ajuste lineal. Mi pregunta es, ¿cómo puedo determinar qué modelo se ajusta mejor a un conjunto de datos particular a partir de la matriz de varianza-covarianza resultante que se devuelve desde el$a * e^{(-b * x)} + c$$a * e^{(-b * x^2)} + c$función scipy.optimize.curve_fit () ? Creo que la variación está en una de las diagonales de esta matriz, pero no estoy seguro de cómo interpretar eso.

ACTUALIZACIÓN: con base en una pregunta similar , espero que la matriz de varianza-covarianza pueda decirme cuál de los tres modelos que estoy intentando ajustar mejor a los datos (estoy tratando de ajustar muchos conjuntos de datos a uno de estos tres modelos).

Las matrices resultantes se ven así para el ejemplo dado:

pcov_lin 
[[  2.02186921e-05  -2.02186920e-04]
 [ -2.02186920e-04   2.76322124e-03]]
pcov_exp
[[  9.05390292e+00  -7.76201283e-02  -9.20475334e+00]
 [ -7.76201283e-02   6.69727245e-04   7.90218415e-02]
 [ -9.20475334e+00   7.90218415e-02   9.36160310e+00]]
pcov_exp_2 
[[  1.38338049e-03  -7.39204594e-07  -7.81208814e-04]
 [ -7.39204594e-07   8.99295434e-09   1.92970700e-06]
 [ -7.81208814e-04   1.92970700e-06   9.14746758e-04]]

Aquí hay un ejemplo de lo que estoy haciendo:

import numpy as np
import matplotlib.pyplot as plt
import scipy as sp
import scipy.optimize

def exp_func(x, a, b, c):
    return a * np.exp(-b * x) + c

def exp_squared_func(x, a, b, c):
    return a * np.exp(-b * x*x*x) + c

def linear_func(x, a, b):
    return a*x + b

def main():
    x = np.array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20], np.float)
    y = np.array([1, 1, 1, 1, 0.805621, 0.798992, 0.84231, 0.728796, 0.819471, 0.570414, 0.355124, 0.276447, 0.159058, 0.0762189, 0.0167807, 0.0118647, 0.000319948, 0.00118267, 0, 0, 0], np.float)

    p0 = [0.7746042467213462, 0.10347274384077858, -0.016253458007293588]
    popt_lin, pcov_lin      = scipy.optimize.curve_fit(linear_func, x, y)
    popt_exp, pcov_exp      = scipy.optimize.curve_fit(exp_func, x, y)
    popt_exp_2, pcov_exp_2  = scipy.optimize.curve_fit(exp_squared_func, x, y)

    plt.figure()
    plt.plot(x, y, 'ko', label="Original data")
    plt.plot(x, linear_func(x, *popt_lin), 'r-', label='linear')
    plt.plot(x, exp_func(x, *popt_exp), 'b-', label='exponential')
    plt.plot(x, exp_squared_func(x, *popt_exp_2), 'g-', label='exponential squared')
    plt.legend()
    plt.show()

if __name__ == '__main__':
    main()

Como aclaración, la variable pcovde scipy.optimize.curve_fites la covarianza estimada de la estimación del parámetro, es decir, en términos generales, dados los datos y un modelo, cuánta información hay en los datos para determinar el valor de un parámetro en el modelo dado. Por lo tanto, realmente no le dice si el modelo elegido es bueno o no. Ver también esto .

El problema de qué es un buen modelo es realmente un problema difícil. Como argumentan los estadísticos

Todos los modelos están equivocados, pero algunos son útiles.

Por lo tanto, el criterio a utilizar en comparación con los diferentes modelos depende de lo que desea lograr.

Por ejemplo, si desea una curva que esté "lo más cerca posible" de los datos, puede seleccionar un modelo que proporcione el residuo más pequeño . En su caso, sería el modelo funcy los parámetros estimados los poptque tienen el valor más bajo al calcular

numpy.linalg.norm(y-func(x, *popt))

Sin embargo, si selecciona un modelo con más parámetros, el residuo disminuirá automáticamente , a costa de una mayor complejidad del modelo. Entonces vuelve a cuál es el objetivo del modelo.