2017-10-12 10 views
1

T 이상 적분 인 함수 f (x)에 데이터를 맞추려고합니다. x는 Integral의 상위 경계선입니다. 나는 scipy.curve_fit()을 사용하여이 작업을 수행하려하지만 curve_fit에 전달할 수있는 함수로 내 적분을 작성하는 방법을 모르겠습니다. 비슷한 문제를 보았지만 내 문제에 맞는 것을 보지 못했습니다.피팅 데이터에 대한 적분

A와 Ea에 대한 추측 값을 제공 할 수 없으므로 지금은 무엇을 할 수 있는지 전혀 알지 못합니다.

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

class get_Ton: 
def __init__(self): 
    self.data=np.genfromtxt('test3.csv', delimiter=',', skip_header=8) 
def loop(self): 
    def fit_tangent(): 
     tck = interpolate.splrep(self.x, self.y, k=2, s=0) 
     dev_1 = interpolate.splev(self.x, tck, der=1) 

     def integrand(T, A, Ea): 
     return A*np.exp(-Ea/(8.314*T)) 
     def polyn(x, A, Ea): 
     return integrate.quad(integrand, 25, x, args=(A, Ea))[0] 

     vcurve = np.vectorize(polyn) 
     p, e = optimize.curve_fit(vcurve, self.x, self.y, [2000, 150]) 
     xd = np.linspace(50, 70, 100) 

     plt.plot(self.x, self.y, "o") 
     plt.plot(xd, vcurve(xd, *p)) 

    for i in range((list(np.shape(self.data)[1:2]))[0]): 
     if i % 2 == 0: 
     self.temp=self.data[:,i] 
     self.scat=self.data[:,i+1] 
     self.x=[26.192, 26.861000000000001, 27.510000000000002, 28.178000000000001, 28.856000000000002, 29.515000000000001, 30.183, 30.823, 31.5, 32.158999999999999, 32.856000000000002, 33.515000000000001, 34.145000000000003, 34.823, 35.491, 36.168999999999997, 36.837000000000003, 37.533999999999999, 38.164000000000001, 38.832000000000001, 39.481000000000002, 40.158999999999999, 40.826999999999998, 41.496000000000002, 42.164000000000001, 42.832000000000001, 43.500999999999998, 44.188000000000002, 44.837000000000003, 45.505000000000003, 46.173000000000002, 46.832000000000001, 47.500999999999998, 48.188000000000002, 48.828000000000003, 49.496000000000002, 50.173999999999999, 50.813000000000002, 51.481000000000002, 52.112000000000002, 52.808999999999997, 53.439, 54.116, 54.765999999999998, 55.453000000000003, 56.101999999999997, 56.761000000000003, 57.429000000000002, 58.078000000000003, 58.737000000000002, 59.442999999999998, 60.082999999999998, 60.770000000000003, 61.448, 62.125999999999998, 62.756, 63.414999999999999, 64.082999999999998, 64.742000000000004, 65.420000000000002, 66.087999999999994, 66.747, 67.415000000000006] 
     self.y=[1553.5, 1595.0, 1497.8, 1695.5999999999999, 1328.7, 1279.0, 1547.8, 1412.8, 1037.0, 1473.5, 1447.4000000000001, 1532.5999999999999, 1484.2, 1169.5, 1395.2, 1183.5999999999999, 949.01999999999998, 1238.0999999999999, 1225.5999999999999, 924.80999999999995, 1650.5999999999999, 803.96000000000004, 1245.7, 1190.0, 1207.0, 1294.0, 1174.9000000000001, 1229.8, 1260.0, 1129.2, 1142.9000000000001, 987.63999999999999, 1389.5999999999999, 1366.0, 1102.0999999999999, 1325.5, 1258.9000000000001, 1285.7, 1217.5, 871.47000000000003, 820.24000000000001, 1388.7, 1391.0, 1400.3, 2482.5999999999999, 3360.5999999999999, 7013.5, 11560.0, 16525.0, 22538.0, 32556.0, 43878.0, 59093.0, 67977.0, 75949.0, 82316.0, 90213.0, 90294.0, 99928.0, 128240.0, 181280.0, 226380.0, 223260.0] 
     fit_tangent() 
     plt.ylim((-100,1000000)) 
     plt.show() 

def main(): 
    this=get_Ton() 
    this.loop() 

if __name__ == "__main__": 
    main() 
+1

"x 위에 적분되는 함수 f (x)"는 의미가 없습니다. 정수의 값은 x에 의존하지 않고 통합 경계에서만 발생합니다. at [25, 65]). 따라서'polyn'는 x에 의존하지 않고'A'와'Ea'를 매개 변수로 취해야하며'optimize.curve_fit'에 의해 최적화 될 것입니다. – user8153

+0

답장을 보내 주셔서 감사합니다.이 신문에 설명 된대로 Arrhenius 변성 모델 FD = 1-exp (-1/B * Int (k (T) dT))를 말합니다 : https : //www.ncbi.nlm .nih.gov/pmc/articles/PMC4709256/ –

+0

귀하의 의견에 따라 질문을 변경했습니다. 고마워. –

답변

2

여기에는 세 가지 문제가 있습니다. 첫째, 함수 polyn은이 변수가 통합되어 있기 때문에 통합 T의 변수에 의존하지 않습니다. 매개 변수 목록에서 T를 제거하십시오. 따라서 trueydata 연산에 수치 중 하나를 드롭 :

trueydata = vcurve(truexdata, 3, 4) 

둘째, quad 튜플 (integral_value, integral_error)을 반환한다. 정수 값만 반환하려면 [0]을 사용하십시오.

def polyn(x, A, Ea): 
    return integrate.quad(integrand, 25, x, args=(A, Ea))[0] 

세 번째로, 매개 변수 값에 대한 초기 추측을 curve_fit으로 제공하십시오. 그렇지 않으면 적합 할 매개 변수의 수를 결정할 수 없다고보고 할 것입니다. 성공하더라도 맹목적으로 초기 추측을 위해 모든 것을 맹목적으로 사용합니다. 최적화 문제의 이해와 함께 인간이 제공 한 초기 추측은 종종 다 변수 최적화의 성공과 실패의 차이입니다.

popt, pcov = optimize.curve_fit(vcurve, xdata, ydata, [2, 3]) 
+0

매우 도움이되는 제안에 대해 감사드립니다. 그에 맞게 스크립트를 수정하고 확장했습니다. 이제 실제 데이터를 추가했습니다. A와 Ea에 대한 추측 값을 제공 할 수 없기 때문에 지금은 무엇을 할 수 있는지 전혀 알지 못합니다. 이제 다른 오류가 나타납니다. ValueError : 둘 이상의 요소가있는 배열의 진리 값이 모호합니다. a.any() 또는 a.all()을 사용하십시오. –

+0

나는 그것을 풀 수 있다고 생각합니다. 그에 따라 코드를 편집했습니다. –