【问题标题】:Exponential fit - sum of exponentials指数拟合 - 指数之和
【发布时间】:2021-05-25 20:42:31
【问题描述】:

我想拟合以下数据:

x(time) y(data)    
0.75;   19.33    
1;      19.04   
1.25;   17.21   
1.5;    12.98   
1.75;   11.59   
2;  9.26   
2.25;   7.66   
2.5;    6.59    
2.75;   5.68    
3;  5.1    
3.25;   4.36   
3.5;    4.43     
3.75;   3.58    
4;  3.01    
4.25;   3.24     
4.5;    3.58     
4.75;   3.13     
5;  3.88     
5.25;   3.19     
5.5;    3.58     
5.75;   3.64   

使用以下代码:

#read text file
data = pd.read_table('episode_5_prova.txt', sep='\t')
#DataFrame
df = pd.DataFrame(data)
#Define your function
def func(x, a, b, c, d, e):
return a*np.exp(-b*x) + c*np.exp(-d*x) + e

#convert dataframe into numpy array
df0=df['time']
x=df0.as_matrix()
df1=df['bi']
y=df1.as_matrix()

# Using the python scipy function curve_fit with function and input variables
popt, pcov = curve_fit(func, x, y)
a, b, c, d, e= popt
fit = func(x, a, b, c, d, e)

fig, ax = plt.subplots()
ax.plot(x, fit, color='r', lw=3)
ax.plot(x, y,'g.')
observed_values=scipy.array(y)
expected_values=scipy.array(fit)
plt.xlim(0,25)
plt.ylim(0,20)

print(a,b,c,d,e)

print(scipy.stats.chisquare(observed_values, f_exp=expected_values, ddof=3))
plt.show()

我得到以下情节: first fit 但是,出于工作的目的,我需要将参数 b 和 c 修复为: b=0.000431062, d=0.000580525 但我没有得到很好的拟合,如下所示: second fit

有人有什么建议吗? 谢谢

【问题讨论】:

  • 您可以将一些文本文件复制/粘贴到问题的正文中,以便我们查看您输入的示例
  • 这样可以吗?谢谢
  • 您发布的数据显示前两个数据点的 Y 值几乎相同,但我在图表上看不到这一点。请您验证一下吗?
  • 是的,对不起,我切入了情节。无论如何,我解决了这个问题......这是关于时间转换的错误......问题是我获得了标准偏差非常高的拟合参数......我不知道如何获得更准确的拟合参数
  • 请您使用新时间转换的更新数据编辑帖子吗?

标签: curve-fitting exponential


【解决方案1】:

这是this question 的副本。

简而言之,您需要一个特征值求解器、一个最小二乘求解器并按照解决方案中列出的步骤进行操作。

Matlab code 可用,您可以将其用作示例以用其他语言(例如 Python)实现。

clear all;
clc;
% get data
dx = 0.02;
x  = (dx:dx:1.5)';
y  =  5*exp(0.5*x) + 4*exp(-3*x) + 2*exp(-2*x) - 3*exp(0.15*x);

% calculate integrals
iy1 = cumtrapz(x, y);
iy2 = cumtrapz(x, iy1);
iy3 = cumtrapz(x, iy2);
iy4 = cumtrapz(x, iy3);

% get exponentials lambdas
Y = [iy1, iy2, iy3, iy4, x.^3, x.^2, x, ones(size(x))];
A = pinv(Y)*y;
%lambdas =
%  -2.9991
%  -1.9997
%   0.5000
%   0.1500

lambdas = eig([A(1), A(2), A(3), A(4); 1, 0, 0, 0; 0, 1, 0, 0; 0, 0, 1, 0]);
lambdas

% get exponentials multipliers
X = [exp(lambdas(1)*x), exp(lambdas(2)*x), exp(lambdas(3)*x), exp(lambdas(4)*x)];
P = pinv(X)*y;
P
%P =
%   4.0042
%   1.9955
%   4.9998
%  -2.9996

【讨论】:

    猜你喜欢
    • 2016-09-16
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 2016-02-13
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    相关资源
    最近更新 更多