¿Existe un algoritmo numérico para encontrar una pendiente asintótica?

23

Tengo una serie de puntos de datos que espero (aproximadamente) sigan una función que asíntota a una línea en grande . Esencialmente, acerca a cero como , y probablemente se puede decir lo mismo de todas las derivadas , , etc. Pero no sé cuál es la forma funcional de f (x) , incluso si tiene una que pueda describirse en términos de funciones elementales.y ( x ) x f ( x ) y ( x ) - ( a x + b ) x f ( x ) f ( x )(xi,yi)y(x)xf(x)y(x)(ax+b)xf(x)f(x)f(x)

Mi objetivo es obtener la mejor estimación posible de la pendiente asintótica a . El método bruto obvio es seleccionar los últimos puntos de datos y hacer una regresión lineal, pero, por supuesto, esto será inexacto si f(x) no se vuelve "lo suficientemente plana" dentro del rango de x para el que tengo datos. El método menos burdo obvio es suponer que f(x)exp(x) (o alguna otra forma funcional particular) y ajustarse a eso usando todos los datos, pero las funciones simples que he probado como exp(x) o 1x no coinciden con los datos en x inferior xdonde f(x)es largo. ¿Existe un algoritmo conocido para determinar la pendiente asintótica que mejoraría, o que podría proporcionar un valor para la pendiente junto con un intervalo de confianza, dada mi falta de conocimiento de cómo se acercan exactamente los datos a la asíntota?


Este tipo de tarea tiende a aparecer con frecuencia en mi trabajo con varios conjuntos de datos, por lo que me interesan principalmente las soluciones generales, pero a pedido, me relaciono con el conjunto de datos en particular que provocó esta pregunta. Como se describe en los comentarios, el algoritmo Wynn ϵ da un valor que, por lo que puedo decir, está algo apagado. Aquí hay una trama:

Datos asintóticamente lineales

(Parece que hay una ligera curva descendente a valores altos de x, pero el modelo teórico para estos datos predice que debería ser asintóticamente lineal).

David Z
fuente
Esto puede ser demasiado elemental, o demasiado vago, para este sitio, pero pensé que la beta privada es el momento de probar tales cosas.
David Z
No, creo que esta es una gran pregunta. No todo debe ser avanzado y elegante. Las buenas soluciones a problemas simples son importantes.
Colin K
@Dan: ¿reemplazar realmente estaba justificado? exp
JM
Tener expertos tiende a hacerme las cosas más difíciles de leer, pero admito que era lo suficientemente pequeño como para no haberlo hecho.
Dan
Realmente no me importa de ninguna manera, solo pensé que podría aprobar las ediciones porque, bueno, ¿por qué no? Obtienes un par de reputación, lo que sea que valga la pena.
David Z

Respuestas:

13

Es un algoritmo bastante aproximado, pero usaría el siguiente procedimiento para una estimación cruda: si, como usted dice, la supuesta que representa su ya es casi lineal a medida que aumenta, lo que yo ' d do es tomar diferencias , y luego usar un algoritmo de extrapolación como la transformación de Shanks para estimar el límite de las diferencias. Con suerte, el resultado es una buena estimación de esta pendiente asintótica.( x i , y i ) x y i + 1 - y if(x)(xi,yi)xyi+1yixi+1xi


Lo que sigue es una demostración de Mathematica . El algoritmo Wynn es una implementación conveniente de la transformación Shanks, y está integrado como la función (oculta) . Probamos el procedimiento en la funciónϵSequenceLimit[]

4x2+3+2x+e4x+3
xdata = RandomReal[{20, 40}, 25];
ydata = Table[(3 + 13*E^(4*x) + 6*E^(4*x)*x + x^2 + 3*E^(4*x)*x^2 + 
      2*E^(4*x)*x^3)/(E^(4*x)*(3 + x^2)), {x, xdata}];

SequenceLimit[Differences[ydata]/Differences[xdata],
              Method -> {"WynnEpsilon", Degree -> 2}]
1.999998

También podría mostrar cuán simple es el algoritmo:

wynnEpsilon[seq_?VectorQ] := 
 Module[{n = Length[seq], ep, res, v, w}, res = {};
  Do[ep[k] = seq[[k]];
   w = 0;
   Do[v = w; w = ep[j];
    ep[j] = 
     v + (If[Abs[ep[j + 1] - w] > 10^-(Precision[w]), ep[j + 1] - w, 
         10^-(Precision[w])])^-1;, {j, k - 1, 1, -1}];
   res = {res, ep[If[OddQ[k], 1, 2]]};, {k, n}];
  Flatten[res]]

Last[wynnEpsilon[Differences[ydata]/Differences[xdata]]]
1.99966

Esta implementación está adaptada del documento de Weniger .

JM
fuente
Es curioso, pero ¿por qué utilizó la forma original de la función, en lugar de combinar todos los términos?
rcollyer
Fue solo con fines de demostración. :) Incluí la expresión -ed para que sepan cuál debería ser la respuesta esperada. Lo que quería era demostrar que se puede hacer este tipo de análisis sobre alguna función complicada de aspecto ...LATEX
JM
¿Qué tan cerca de plano deben estar los puntos para que el algoritmo sea efectivo?
rcollyer
2
Bien, última pregunta (lo juro), ¿puedes generar un error vinculado a la estimación?
rcollyer
1
Eso es un poco más complicado. He visto algunos métodos sugeridos en algunos documentos, pero confieso que no he hecho experimentos con ellos. (Tal vez debería hacerlo, uno de estos días.) El libro de Brezinski y Redivo-Zaglia tiene algunos consejos que tal vez quieras analizar.
JM