2009-06-22 4 views
11

Есть ли лучший способ найти который X дает мне Y Я ищу в SciPy? Я только начал использовать SciPy, и я не очень хорошо знаком с каждой функцией.Интерполяция в SciPy: Обнаружение X, который производит Y

import numpy as np 
import matplotlib.pyplot as plt 
from scipy import interpolate 

x = [70, 80, 90, 100, 110] 
y = [49.7, 80.6, 122.5, 153.8, 163.0] 
tck = interpolate.splrep(x,y,s=0) 
xnew = np.arange(70,111,1) 
ynew = interpolate.splev(xnew,tck,der=0) 
plt.plot(x,y,'x',xnew,ynew) 
plt.show() 
t,c,k=tck 
yToFind = 140 
print interpolate.sproot((t,c-yToFind,k)) #Lowers the spline at the abscissa 
+0

Вы можете уточнить, что вы хотите быть лучше? Производительность, точность, краткость? – Ryan

ответ

16

Класс UnivariateSpline в scipy делает выполнение сплайнов намного более питоническим.

x = [70, 80, 90, 100, 110] 
y = [49.7, 80.6, 122.5, 153.8, 163.0] 
f = interpolate.UnivariateSpline(x, y, s=0) 
xnew = np.arange(70,111,1) 

plt.plot(x,y,'x',xnew,f(xnew)) 

Чтобы найти х в точке у, то сделать:

yToFind = 140 
yreduced = np.array(y) - yToFind 
freduced = interpolate.UnivariateSpline(x, yreduced, s=0) 
freduced.roots() 

Я думал интерполирования х через у может работать, но это занимает несколько иной маршрут. Это может быть ближе к большему количеству очков.

+1

Разве это не потребует в два раза больше вычислений ЦП, поскольку они интерполируют практически тот же набор данных два раза? – JcMaco

+0

@JcMaco, первое использование UnivariateSpline - просто сделать красивый сюжет. Второе использование - это то, что фактически дает значения. – Theran

+0

Крейг прав, можете ли вы исправить его на своем примере, так как это здорово! –

1

Если все, что вам нужно, это линейная интерполяция, вы можете использовать функцию interp в NumPy.

+0

Я бы предпочел сплайн-интерполяцию. Как бы функция интерпретации помогла мне легче решить мою проблему? – JcMaco

+0

В вашем вопросе не указано, какой тип интерполяции вам нужен - если линейка недостаточно хороша для вашей проблемы, тогда я не думаю, что интерпретация поможет. –

0

Возможно, я неправильно понял ваш вопрос, если так мне жаль. Я не думаю, что вам нужно использовать SciPy. NumPy имеет функцию наименьших квадратов.

#!/usr/bin/env python 

from numpy.linalg.linalg import lstsq 



def find_coefficients(data, exponents): 
    X = tuple((tuple((pow(x,p) for p in exponents)) for (x,y) in data)) 
    y = tuple(((y) for (x,y) in data)) 
    x, resids, rank, s = lstsq(X,y) 
    return x 

if __name__ == "__main__": 
    data = tuple((
     (1.47, 52.21), 
     (1.50, 53.12), 
     (1.52, 54.48), 
     (1.55, 55.84), 
     (1.57, 57.20), 
     (1.60, 58.57), 
     (1.63, 59.93), 
     (1.65, 61.29), 
     (1.68, 63.11), 
     (1.70, 64.47), 
     (1.73, 66.28), 
     (1.75, 68.10), 
     (1.78, 69.92), 
     (1.80, 72.19), 
     (1.83, 74.46) 
    )) 
    print find_coefficients(data, range(3)) 

Это вернет [128.81280358 -143.16202286 61.96032544].

>>> x=1.47 # the first of the input data 
>>> 128.81280358 + -143.16202286*x + 61.96032544*(x**2) 
52.254697219095988 

0,04 из, не плохо

+0

Проблема заключается в том, какой номер даст мне 52.21. Конечно, может быть много решений, если интерполяция квадратична (или более высокая мощность). – JcMaco