scipy에서 제공하는 보간기를 사용하기 위해 수학자 동료가 개발 한 수동 보간기를 사용하는 프로그램을 이식하려고합니다. scipy 보간기를 사용하거나 래핑하여 가능한 한 이전 보간기에 가깝게 작동하도록하고 싶습니다.
두 기능의 주요 차이점은 원래 보간 기에서 입력 값이 입력 범위보다 크거나 작 으면 원래 보간 기가 결과를 추정한다는 것입니다. scipy 보간기로 이것을 시도하면 ValueError. 이 프로그램을 예로 고려하십시오.
import numpy as np
from scipy import interpolate
x = np.arange(0,10)
y = np.exp(-x/3.0)
f = interpolate.interp1d(x, y)
print f(9)
print f(11) # Causes ValueError, because it's greater than max(x)
충돌하는 대신 최종 라인이 단순히 선형 외삽을 수행하여 첫 번째와 마지막 두 점으로 정의 된 기울기를 무한대로 계속하도록 만드는 합리적인 방법이 있습니까?
실제 소프트웨어에서는 실제로 exp 함수를 사용하지 않습니다. 여기에서는 설명 용으로 만 사용합니다!
답변
1. 상수 외삽
interpscipy의 함수를 사용할 수 있으며 왼쪽 및 오른쪽 값을 범위를 초과하는 상수로 외삽합니다.
>>> from scipy import interp, arange, exp
>>> x = arange(0,10)
>>> y = exp(-x/3.0)
>>> interp([9,10], x, y)
array([ 0.04978707, 0.04978707])
2. 선형 (또는 기타 사용자 정의) 외삽
선형 외삽을 처리하는 보간 함수 주위에 래퍼를 작성할 수 있습니다. 예를 들면 :
from scipy.interpolate import interp1d
from scipy import arange, array, exp
def extrap1d(interpolator):
xs = interpolator.x
ys = interpolator.y
def pointwise(x):
if x < xs[0]:
return ys[0]+(x-xs[0])*(ys[1]-ys[0])/(xs[1]-xs[0])
elif x > xs[-1]:
return ys[-1]+(x-xs[-1])*(ys[-1]-ys[-2])/(xs[-1]-xs[-2])
else:
return interpolator(x)
def ufunclike(xs):
return array(list(map(pointwise, array(xs))))
return ufunclike
extrap1d보간 함수를 받고 외삽 할 수있는 함수를 반환합니다. 다음과 같이 사용할 수 있습니다.
x = arange(0,10)
y = exp(-x/3.0)
f_i = interp1d(x, y)
f_x = extrap1d(f_i)
print f_x([9,10])
산출:
[ 0.04978707 0.03009069]
답변
InterpolatedUnivariateSpline을 살펴볼 수 있습니다.
다음은 그것을 사용하는 예입니다.
import matplotlib.pyplot as plt
import numpy as np
from scipy.interpolate import InterpolatedUnivariateSpline
# given values
xi = np.array([0.2, 0.5, 0.7, 0.9])
yi = np.array([0.3, -0.1, 0.2, 0.1])
# positions to inter/extrapolate
x = np.linspace(0, 1, 50)
# spline order: 1 linear, 2 quadratic, 3 cubic ...
order = 1
# do inter/extrapolation
s = InterpolatedUnivariateSpline(xi, yi, k=order)
y = s(x)
# example showing the interpolation for linear, quadratic and cubic interpolation
plt.figure()
plt.plot(xi, yi)
for order in range(1, 4):
s = InterpolatedUnivariateSpline(xi, yi, k=order)
y = s(x)
plt.plot(x, y)
plt.show()
답변
SciPy 버전 0.17.0부터 외삽을 허용 하는 scipy.interpolate.interp1d에 대한 새로운 옵션 이 있습니다. 호출에서 fill_value = ‘extrapolate’를 설정하기 만하면됩니다. 이 방법으로 코드를 수정하면 다음이 제공됩니다.
import numpy as np
from scipy import interpolate
x = np.arange(0,10)
y = np.exp(-x/3.0)
f = interpolate.interp1d(x, y, fill_value='extrapolate')
print f(9)
print f(11)
출력은 다음과 같습니다.
0.0497870683679
0.010394302658
답변
scipy.interpolate.splrep은 어떻습니까 (차수가 1이고 평활화 없음) :
>> tck = scipy.interpolate.splrep([1, 2, 3, 4, 5], [1, 4, 9, 16, 25], k=1, s=0)
>> scipy.interpolate.splev(6, tck)
34.0
34 = 25 + (25-16)이므로 원하는 것을 수행하는 것 같습니다.
답변
다음은 numpy 패키지 만 사용하는 대체 방법입니다. numpy의 배열 함수를 활용하므로 큰 배열을 보간 / 외삽 할 때 더 빠를 수 있습니다.
import numpy as np
def extrap(x, xp, yp):
"""np.interp function with linear extrapolation"""
y = np.interp(x, xp, yp)
y = np.where(x<xp[0], yp[0]+(x-xp[0])*(yp[0]-yp[1])/(xp[0]-xp[1]), y)
y = np.where(x>xp[-1], yp[-1]+(x-xp[-1])*(yp[-1]-yp[-2])/(xp[-1]-xp[-2]), y)
return y
x = np.arange(0,10)
y = np.exp(-x/3.0)
xtest = np.array((8.5,9.5))
print np.exp(-xtest/3.0)
print np.interp(xtest, x, y)
print extrap(xtest, x, y)
편집 : Mark Mikofski가 제안한 “extrap”기능 수정 :
def extrap(x, xp, yp):
"""np.interp function with linear extrapolation"""
y = np.interp(x, xp, yp)
y[x < xp[0]] = yp[0] + (x[x<xp[0]]-xp[0]) * (yp[0]-yp[1]) / (xp[0]-xp[1])
y[x > xp[-1]]= yp[-1] + (x[x>xp[-1]]-xp[-1])*(yp[-1]-yp[-2])/(xp[-1]-xp[-2])
return y
답변
그것은 사용하기 빨라질 수 부울 인덱싱 과 대용량 데이터 부울 인덱싱 쉽고 빠르게 비교를 허용하는 반면 모든 점은, 간격 외측에있는 경우, 알고리즘 검사 이후.
예를 들면 :
# Necessary modules
import numpy as np
from scipy.interpolate import interp1d
# Original data
x = np.arange(0,10)
y = np.exp(-x/3.0)
# Interpolator class
f = interp1d(x, y)
# Output range (quite large)
xo = np.arange(0, 10, 0.001)
# Boolean indexing approach
# Generate an empty output array for "y" values
yo = np.empty_like(xo)
# Values lower than the minimum "x" are extrapolated at the same time
low = xo < f.x[0]
yo[low] = f.y[0] + (xo[low]-f.x[0])*(f.y[1]-f.y[0])/(f.x[1]-f.x[0])
# Values higher than the maximum "x" are extrapolated at same time
high = xo > f.x[-1]
yo[high] = f.y[-1] + (xo[high]-f.x[-1])*(f.y[-1]-f.y[-2])/(f.x[-1]-f.x[-2])
# Values inside the interpolation range are interpolated directly
inside = np.logical_and(xo >= f.x[0], xo <= f.x[-1])
yo[inside] = f(xo[inside])
제 경우에는 300000 포인트의 데이터 세트를 사용하면 25.8 초에서 0.094 초로 속도가 빨라지고 250 배 이상 빠릅니다 .
답변
초기 배열에 포인트를 추가하여 수행했습니다. 이런 식으로 나는 스스로 만든 함수를 정의하는 것을 피하고 선형 외삽 (아래 예제에서 : 오른쪽 외삽)은 괜찮아 보입니다.
import numpy as np
from scipy import interp as itp
xnew = np.linspace(0,1,51)
x1=xold[-2]
x2=xold[-1]
y1=yold[-2]
y2=yold[-1]
right_val=y1+(xnew[-1]-x1)*(y2-y1)/(x2-x1)
x=np.append(xold,xnew[-1])
y=np.append(yold,right_val)
f = itp(xnew,x,y)
