2011-09-08 429 views
13

Soy bastante nuevo en la programación y pensé en intentar escribir una función de interpolación lineal.Interpolación lineal - Python

decir que estoy datos dados como sigue: x = [1, 2,5, 3,4, 5,8, 6] y = [2, 4, 5,8, 4,3, 4]

Quiero diseñar una función que se interpolará linealmente entre 1 y 2.5, 2.5 a 3.4 y así sucesivamente utilizando Python.

He intentado mirar a través de http://docs.python.org/tutorial/ pero sigo sin poder entenderlo.

+0

Esto es ... no es fácil. ¿Qué has intentado? – zellio

+0

-1 como demasiado general. usted no entiende cómo programar, o cómo hacer el algoritmo en Python? – steabert

+0

Bueno, siendo un nuevo aprendiz, me he lanzado a lo más profundo, por así decirlo. Estaba pensando en usar sentencias 'for' o 'if' en un algoritmo. Entonces entre numerosos rangos de x. – Helpless

Respuesta

13

Según entiendo su pregunta, ¿desea escribir alguna función y = interpolate(x_values, y_values, x), que le dará el valor y en unos x? La idea básica luego sigue estos pasos:

  1. Encuentra los índices de los valores en x_values que definen un intervalo que contiene x. Por ejemplo, para x=3 con sus listas de ejemplo, el intervalo que contiene sería [x1,x2]=[2.5,3.4], y los índices sería i1=1, i2=2
  2. calcular la pendiente en este intervalo por (y_values[i2]-y_values[i1])/(x_values[i2]-x_values[i1]) (es decir dy/dx).
  3. El valor en x ahora es el valor en x1 más la pendiente multiplicada por la distancia desde x1.

Usted, además, tendrá que decidir qué ocurre si x está fuera del intervalo de x_values, o bien se trata de un error, o se puede interpolar "hacia atrás", suponiendo que la pendiente es la misma que la primera último intervalo /.

¿Te ayudó o necesitaste más consejos específicos?

+0

No, eso es perfecto, ¡muchas gracias! – Helpless

+3

No puede ser correcto restar y_values ​​[i2] de sí mismo. ¿Debería ser '(y_values ​​[i2] -y_values ​​[i1])'? –

+3

@MartinBurch: Unos años más tarde, pero ... ¡Gracias, corregido! :) – carlpett

10

me ocurrió una solución bastante elegante (en mi humilde opinión), por lo que no puede resistirse a publicarla:

from bisect import bisect_left 

class Interpolate(object): 
    def __init__(self, x_list, y_list): 
     if any(y - x <= 0 for x, y in zip(x_list, x_list[1:])): 
      raise ValueError("x_list must be in strictly ascending order!") 
     x_list = self.x_list = map(float, x_list) 
     y_list = self.y_list = map(float, y_list) 
     intervals = zip(x_list, x_list[1:], y_list, y_list[1:]) 
     self.slopes = [(y2 - y1)/(x2 - x1) for x1, x2, y1, y2 in intervals] 

    def __getitem__(self, x): 
     i = bisect_left(self.x_list, x) - 1 
     return self.y_list[i] + self.slopes[i] * (x - self.x_list[i]) 

I mapa para float de manera que la división entera (pitón < = 2,7) no entrará en funcionamiento y arruinar las cosas si x1, x2, y1 y y2 son enteros para algunos iterval.

En __getitem__ estoy tomando ventaja del hecho de que self.x_list se ordena en orden ascendente utilizando bisect_left a (muy) encontrar rápidamente el índice del elemento más grande más pequeño que x en self.x_list.

Uso de la clase como esta:

i = Interpolate([1, 2.5, 3.4, 5.8, 6], [2, 4, 5.8, 4.3, 4]) 
# Get the interpolated value at x = 4: 
y = i[4] 

yo no he tratado con las condiciones de borde a todos aquí, por razones de simplicidad. Tal como está, i[x] para x < 1 funcionará como si la línea desde (2.5, 4) a (1, 2) se hubiera extendido hasta menos el infinito, mientras que i[x] para x == 1 o x > 6 levantará un IndexError. Sería mejor plantear un IndexError en todos los casos, pero esto se deja como un ejercicio para el lector.:)

+1

Me gustaría usar '' __call__' en lugar de '__getitem__' para ser preferible en general, generalmente es una función de interpolación *. – Dave

23
import scipy.interpolate 
y_interp = scipy.interpolate.interp1d(x, y) 
print y_interp(5.0) 

scipy.interpolate.interp1d hace la interpolación lineal por y se puede personalizar para manejar las condiciones de error.

1

Su solución no funcionó en Python 2.7. Hubo un error al verificar el orden de los elementos x. Tuve que cambiar a este código para conseguir que funcione:

from bisect import bisect_left 
class Interpolate(object): 
    def __init__(self, x_list, y_list): 
     if any([y - x <= 0 for x, y in zip(x_list, x_list[1:])]): 
      raise ValueError("x_list must be in strictly ascending order!") 
     x_list = self.x_list = map(float, x_list) 
     y_list = self.y_list = map(float, y_list) 
     intervals = zip(x_list, x_list[1:], y_list, y_list[1:]) 
     self.slopes = [(y2 - y1)/(x2 - x1) for x1, x2, y1, y2 in intervals] 
    def __getitem__(self, x): 
     i = bisect_left(self.x_list, x) - 1 
     return self.y_list[i] + self.slopes[i] * (x - self.x_list[i]) 
1

En lugar de la extrapolación de los extremos, se podría volver la extensión de la y_list. La mayoría de las veces su aplicación se comporta bien, y el Interpolate[x] estará en el x_list. Los (presumiblemente) efectos lineales de la extrapolación de los extremos pueden inducir a error a creer que sus datos se comportan bien.

  • Volviendo resultado no lineal (delimitada por el contenido de x_list y y_list) el comportamiento de su programa que puede alertar a un problema de valores en gran medida fuera x_list. (Comportamiento lineal va los plátanos cuando se les da entradas no lineales!)

  • Volviendo la extensión de la y_list para Interpolate[x] fuera del x_list también significa que se conoce el alcance de su valor de salida. Si extrapola en base a x mucho, mucho menos que x_list[0] o x mucho, mucho más que x_list[-1], su resultado de devolución podría estar fuera del rango de valores que esperaba.

    def __getitem__(self, x): 
        if x <= self.x_list[0]: 
         return self.y_list[0] 
        elif x >= self.x_list[-1]: 
         return self.y_list[-1] 
        else: 
         i = bisect_left(self.x_list, x) - 1 
         return self.y_list[i] + self.slopes[i] * (x - self.x_list[i]) 
    
+0

Me gustaría usar '' __call__' en lugar de '__getitem__' para ser preferible en general, por lo general es una función de interpolación *. – Dave

Cuestiones relacionadas