I need to fits the attached data as in image

3 次查看（过去 30 天）
Mushtaq Al-Jubbori
评论： Mushtaq Al-Jubbori 2021-6-10
I need to fit the attached data as in the image below:
  2 个评论
Scott MacKenzie
Scott MacKenzie 2021-6-5
Have you made any attempt yourself to do this? If no, then perhaps you should. If yes, then it might help if you showed the code you've put together so far.

请先登录，再进行评论。

请先登录，再回答此问题。

采纳的回答

Alex Sha
Alex Sha 2021-6-6
Hi, the fitting function "y = p1/(1+Exp(p2+p3*x)+p4*x)^p5" is pretty good for all data set, where y=L, x=T
1: T1-L1:
Root of Mean Square Error (RMSE): 0.0481024350848356
Sum of Squared Residual: 0.0277661311330899
Correlation Coef. (R): 0.999085841980414
R-Square: 0.998172519645713
Parameter Best Estimate
---------- -------------
p1 6.26675447680147
p2 109.105635295479
p3 -67.1792841156261
p4 404935160.81852
p5 0.0161845071634169
2: T2-L2:
Root of Mean Square Error (RMSE): 0.0728301819569392
Sum of Squared Residual: 0.111388943481498
Correlation Coef. (R): 0.999334306206347
R-Square: 0.998669055560921
Parameter Best Estimate
---------- -------------
p1 7.42506477062961
p2 19.1556613661477
p3 -9.26234046573527
p4 -0.0219401101448161
p5 0.101288835010354
3: T3-L3:
Root of Mean Square Error (RMSE): 0.134092463526771
Sum of Squared Residual: 0.413558141817605
Correlation Coef. (R): 0.999276494987953
R-Square: 0.998553513435409
Parameter Best Estimate
---------- -------------
p1 11.4516493147292
p2 17.1095166703901
p3 -4.66159981573592
p4 0.0192821696183541
p5 0.129232672490336
4: T4-L4:
Root of Mean Square Error (RMSE): 0.174646688455142
Sum of Squared Residual: 1.06755130259216
Correlation Coef. (R): 0.999289186413636
R-Square: 0.998578878083226
Parameter Best Estimate
---------- -------------
p1 15.0727185425809
p2 61.2119358895463
p3 -10.7637556382665
p4 0.236072827848106
p5 0.037477078210382
5: T5-L5:
Root of Mean Square Error (RMSE): 0.154312409495508
Sum of Squared Residual: 0.500058714210495
Correlation Coef. (R): 0.999654251549009
R-Square: 0.999308622640009
Parameter Best Estimate
---------- -------------
p1 17.5989211228526
p2 138.340503536049
p3 -17.2818282854456
p4 -0.0357446437398955
p5 0.0179266207350469
  11 个评论
显示 9更早的评论
Image Analyst
Image Analyst 2021-6-10
@Mushtaq Al-Jubbori, you say "I wish to developement the function coth(x)-1/x"
Well, I gave you code for that here:
https://www.mathworks.com/matlabcentral/answers/848200-i-need-to-fits-the-attached-data-as-in-image#comment_1569470
You're free to develop the function further if you want.
Mushtaq Al-Jubbori
Mushtaq Al-Jubbori 2021-6-10
Thank you, but when run the cod appear.. Undefined function or variable 'FitCoth'

请先登录，再进行评论。

更多回答（5 个）

Sulaymon Eshkabilov
You can start using curve fitting toolbox, cftool that is quite straightforward and does not require any addional coding.
An alternative way is the least squares method for that you have you data ready. You'd need to generate Vandermonde matrix and then compute the fit model coefficients.
  4 个评论
显示 2更早的评论
John D'Errico
John D'Errico 2021-6-6
NO. You do NOT want to use a polynomial model for that data!
Polynomials NEVER have the property that they roll over and then become flat.
Sulaymon Eshkabilov
Sulaymon Eshkabilov 2021-6-8
@John D'Errico This polynomial model is shown here mean to be an example not the proposed model. For V matrix it is viable to insert exp(), sin(), cos() tan() and polynonials in combination.

请先登录，再进行评论。

Image Analyst
Image Analyst 2021-6-6
@Mushtaq Al-Jubbori. What I did was to find where the flat part starts and fit the data to the left of that to an exponential growth curve. As an output you have the location where the flat part starts and the coefficients of the exponential fitted equation. Try this (click the copy icon in the upper right of the code block below and then paste into MATLAB):
clc; % Clear command window.
fprintf('Running %s.m ...\n', mfilename);
clear; % Delete all variables.
close all; % Close all figure windows except those created by imtool.
workspace; % Make sure the workspace panel is showing.
fontSize = 14;
L1=[1.4 1.9 2.4 3.132 4.2 4.532 4.532 4.4 4.532 4.532 4.4 4.5];%1.5Mev
L2=[1.65 2.2 2.8 3.4 4.266 5.6 6.6 7.4 7.5 7.534 7.5 7.5 7.4 7.4 7.4 7.5 7.6 7.532 7.4 7.6 7.6];%2.26MeV
L3=[1.8 2.2 2.8 3 3.8 4.8 5.8 6.4 7.8 8.6 9.8 10.6 11.2 11.4 11.2 11.4 11.2 11.4 11 11.2 11.4 11.2 11.4 ];%3.2MeV
L4=[2.05 2.4 2.8 3.4 3.7 4 4.5 5 6.3 7.2 8 8.6 9.8 10.4 11.2 12.2 13.8 14.6 14.6 14.532 14.532 14.6 14.6 14.532 14.532 14.4 14.4 14.4 14.4 14.532 14.5 14.2 14.4 14.4 14.6];%4MeV
L5=[3 3.4 3.8 4.4 5.2 6 6.8 8 9.2 10.8 12.8 15 16.4 17.4 17.8 17.6 17.8 17.8 17.6 17.8 17.8];%4.85MeV
T1=0.25:0.25:3;
T2=[0.5:0.25:2.25 2.75:0.25:5 5.5 5.75 6];
T3=[0.75:0.25:1.75 2.25:0.25:4 4.5:0.25:5 5.5 6 6.25 6.75:0.25:7.5 ];
T4=[1 1.25 1.75:0.25:3 3.5:0.25:9.25 9.75 10.25 10.75];
T5=[2:0.5:7.5 7.7 8:0.5:10.5 11.5 12];
subplot(5, 2, 1);
plot(T1,L1)
grid on;
title('L1 vs. T1', 'FontSize', fontSize);
kneeIndex = FindKneeIndex(L1, T1);
xline(T1(kneeIndex), 'Color', 'r', 'LineWidth', 2);
% Now fit the left portion of the data that follows an exponential to an exponential formula.
coefficients1 = FitExponentialGrowth(T1(1:kneeIndex), L1(1:kneeIndex), 2)
subplot(5, 2, 3);
plot(T2,L2)
grid on;
title('L2 vs. T2', 'FontSize', fontSize);
kneeIndex = FindKneeIndex(L2, T2);
xline(T2(kneeIndex), 'Color', 'r', 'LineWidth', 2);
% Now fit the left portion of the data that follows an exponential to an exponential formula.
coefficients2 = FitExponentialGrowth(T2(1:kneeIndex), L2(1:kneeIndex), 4)
subplot(5, 2, 5);
plot(T3,L3)
grid on;
title('L3 vs. T3', 'FontSize', fontSize);
kneeIndex = FindKneeIndex(L3, T3);
xline(T3(kneeIndex), 'Color', 'r', 'LineWidth', 2);
% Now fit the left portion of the data that follows an exponential to an exponential formula.
coefficients3 = FitExponentialGrowth(T3(1:kneeIndex), L3(1:kneeIndex), 6)
subplot(5, 2, 7);
plot(T4,L4)
grid on;
title('L4 vs. T4', 'FontSize', fontSize);
kneeIndex = FindKneeIndex(L4, T4);
xline(T4(kneeIndex), 'Color', 'r', 'LineWidth', 2);
% Now fit the left portion of the data that follows an exponential to an exponential formula.
coefficients4 = FitExponentialGrowth(T4(1:kneeIndex), L4(1:kneeIndex), 8)
subplot(5, 2, 9);
plot(T5,L5)
grid on;
title('L5 vs. T5', 'FontSize', fontSize);
kneeIndex = FindKneeIndex(L5, T5);
xline(T5(kneeIndex), 'Color', 'r', 'LineWidth', 2);
% Now fit the left portion of the data that follows an exponential to an exponential formula.
coefficients5 = FitExponentialGrowth(T5(1:kneeIndex), L5(1:kneeIndex), 10)
% Maximize figure window.
g = gcf;
g.WindowState = 'maximized';
fprintf('Done running %s.m\n', mfilename);
function kneeIndex = FindKneeIndex(L, T)
slopes = (L(end) - L) ./ (T(end) - T);
kneeIndex = find(slopes < 0.15, 1, 'first');
% figure;
% plot(T, slopes, 'b.-');
% grid on;
end
function coefficients = FitExponentialGrowth(X, Y, plotNumber)
fontSize = 14;
% Convert X and Y into a table, which is the form fitnlm() likes the input data to be in.
tbl = table(X(:), Y(:));
% Define the model as Y = a + exp(-b*x)
% Note how this "x" of modelfun is related to big X and big Y.
% x((:, 1) is actually X and x(:, 2) is actually Y - the first and second columns of the table.
modelfun = @(b,x) b(1) * exp(b(2)*x(:, 1)) + b(3);
beta0 = [11, .5, 0]; % Guess values to start with. Just make your best guess.
% Now the next line is where the actual model computation is done.
mdl = fitnlm(tbl, modelfun, beta0);
% Now the model creation is done and the coefficients have been determined.
% YAY!!!!
% Extract the coefficient values from the the model object.
% The actual coefficients are in the "Estimate" column of the "Coefficients" table that's part of the mode.
coefficients = mdl.Coefficients{:, 'Estimate'}
% Create smoothed/regressed data using the model:
yFitted = coefficients(1) * exp(coefficients(2)*X) + coefficients(3);
% Now we're done and we can plot the smooth model as a red line going through the noisy blue markers.
subplot(5, 2, plotNumber);
plot(X, Y, 'b.', 'MarkerSize', 20);
hold on;
plot(X, yFitted, 'r-', 'LineWidth', 2);
grid on;
title('Exponential Regression with fitnlm()', 'FontSize', fontSize);
xlabel('X', 'FontSize', fontSize);
ylabel('Y', 'FontSize', fontSize);
legendHandle = legend('Noisy Y', 'Fitted Y', 'Location', 'northwest');
legendHandle.FontSize = 12;
formulaString = sprintf('Y = %.3f * exp(%.3f * X) + %.3f', coefficients(1), coefficients(2), coefficients(3))
text(5, 17000, formulaString, 'FontSize', 12, 'FontWeight', 'bold');
end
  12 个评论
显示 10更早的评论
Image Analyst
Image Analyst 2021-6-8
编辑：Image Analyst 2021-6-8
I don't know if coth(x)-1/x is ok or not. It doesn't seem to fit as well as the piecewise exponential growth function I suggested.
x = linspace(-100, 100, 1000);
y = coth(x) - 1 ./ x;
plot(x, y, 'b-', 'lineWidth', 2);
grid on;
but let's say that YOU say it's ok. So please explain the physical meaning of the curve, meaning how it theoretically applies to the physical process you are modeling, which is "Track length of alpha pariclese and etchin time". Presumably someone published a paper where the physics says it should be that equation theoretically.
I don't think the fit looks very good. Defined the model as Y = a * coth(b * x) - c ./ x
Code is attached.
Mushtaq Al-Jubbori
Mushtaq Al-Jubbori 2021-6-8
Did we can develobment the function by + or - any parameters?

请先登录，再进行评论。

Scott MacKenzie
Scott MacKenzie 2021-6-6
编辑：Scott MacKenzie 2021-6-6
Here's what I put together. This is only for the T5-L5 data. You can repeat this for the rest of the data and build up the plot with the remaining points and lines.
L5=[3 3.4 3.8 4.4 5.2 6 6.8 8 9.2 10.8 12.8 15 16.4 17.4 17.8 17.6 17.8 17.8 17.6 17.8 17.8];%4.85MeV
T5=[2:0.5:7.5 7.7 8:0.5:10.5 11.5 12];
threshold = 0.2; % adjust accordingly
% find elbow (index where curve transitions to flat line)
elbow = find(diff(L5) < threshold, 1);
% plot markers, set axis limits
plot(T5, L5, 'or', 'markerfacecolor', 'r');
ax = gca;
ax.YLim = [0 20];
ax.XLim = [0 15];
% plot flat line
x1 = T5(elbow); x2 = T5(end);
y1 = mean(L5(elbow):L5(end)); y2 = y1;
line([x1 x2], [y1 y2], 'color', 'r', 'linewidth', 1.5);
hold on;
% plot curved line using polynomial fit
x = T5(1:elbow);
ySample = L5(1:elbow);
pf = polyfit(x,ySample,5); % order 5 looks pretty good
y = polyval(pf, x)
y(end) = L5(elbow); % ensure curved line meets flat line
plot(x, y, 'r', 'linewidth', 1.5);
Sulaymon Eshkabilov
This nonlinear least squares fit gives quite nice fit models:
...
x = T2;
y = L2;
p = [1, 3]; % Guess: Fit model parameters
q = [-1, -3]; % Guess: Fit model parameters
plot(x,y,'r*')
xlabel 't'
ylabel 'Response'
p = optimvar('p',2);
q = optimvar('q',2);
fun = p(1)*(1+exp(q(1)+q(2)*x) + p(2)*x).^-0.65;
obj = sum((fun - y).^2);
LSQ_Prob = optimproblem("Objective",obj);
x0.p = [1/2, 5/2];
x0.q = [-1/2,-5/2];
show(LSQ_Prob) % Display the problem formulation and parameters
[Sol,fval] = solve(LSQ_Prob,x0) % Show the results
figure
FM_val = evaluate(fun,Sol);
plot(x,y,'r*', x, FM_val,'b-')
legend('Original Data','Fitted Curve', 'location', 'southeast')
xlabel('T')
ylabel('L & Fit Model')
title("Data vs. Fitted Model Response")
SStot = sum((y-mean(y)).^2);
SSres = sum((y-FM_val).^2);
Rsq = 1 - SSres/SStot;
gtext(['R^2 = ' num2str(Rsq)])
gtext(['Fit model: ' num2str(Sol.p(1)) '*(1+exp(' num2str(Sol.q(1)) ' + ' num2str(Sol.q(2)) '*x) + ' num2str(Sol.p(2)) '*x).^{-0.65}'])
fprintf('Found fit model is: %f *(1+exp(%f+%f*x) + %f*x)^-0.65 \n', [Sol.p(1), Sol.q(1), Sol.q(2), Sol.p(2)])
  7 个评论
显示 5更早的评论
Mushtaq Al-Jubbori
Mushtaq Al-Jubbori 2021-6-8
Very thanks but function coth not ok, there for I cant not explin else the function is complet (no. Parameters)
Image Analyst
Image Analyst 2021-6-8
OK, fine. Then use the piecewise function I gave you, or use ad hoc function mysteriously found by that third party software.
https://www.mathworks.com/matlabcentral/answers/848200-i-need-to-fits-the-attached-data-as-in-image#answer_718230
but you may need to upgrade your MATLAB, as you said.

请先登录，再进行评论。

Mushtaq Al-Jubbori
Thanks
Sulaymon Eshkabilov
But the cod dont run in Matlab-2016, the function optimvar undefine

类别

在 Help Center 和 File Exchange 中查找有关 Model Building and Assessment 的更多信息

标签

产品

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by