I need to fits the attached data as in image

I need to fit the attached data as in the image below:

[Scott MacKenzie]

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.

[Mushtaq Al-Jubbori]

Yes, I used tanh(x) but I need another function?

Answer 1

[Alex Sha]

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

[Image Analyst]

Looks like that model fits pretty well. Where does that model come from? Is it the theoretical model for some actual real world physical process?

[Scott MacKenzie]

@Alex Sha Could you please post the code you used to generate one of the above fits and plots.

[Alex Sha]

Hi, the function was got automatically through another package of 1stOpt. Actually, the original function of "y=A1*tanh(U)^A4, where U=exp((x-A2)/A3)" is also good enough:

1: T1-L1:

Root of Mean Square Error (RMSE): 0.0597205155644429
Sum of Squared Residual: 0.0427984797513945
Correlation Coef. (R): 0.998591934579049
R-Square: 0.997185851806327

Parameter	Best Estimate
----------	-------------
a1	4.49354483105404
a2	1.2952779431858
a3	0.182566357732861
a4	0.204961211278745

2: T2-L2:

Root of Mean Square Error (RMSE): 0.0787759184709166
Sum of Squared Residual: 0.130318551949667
Correlation Coef. (R): 0.999222163002948
R-Square: 0.998444931036289

Parameter	Best Estimate
----------	-------------
a1	7.50386733656648
a2	2.02357393772506
a3	0.425853364237124
a4	0.421653184517995

3: T3-L3:

Root of Mean Square Error (RMSE): 0.130459406471193
Sum of Squared Residual: 0.391452104946764
Correlation Coef. (R): 0.99931463290777
R-Square: 0.99862973554359

Parameter	Best Estimate
----------	-------------
a1	11.2849810802541
a2	3.62185746086777
a3	0.583978847866495
a4	0.356432637089129

4: T4-L4:

Root of Mean Square Error (RMSE): 0.17944782733443
Sum of Squared Residual: 1.12705329572666
Correlation Coef. (R): 0.999239734136879
R-Square: 0.998480046277941

Parameter	Best Estimate
----------	-------------
a1	14.4986426026826
a2	5.50830739927991
a3	0.561356452592117
a4	0.235099297830832

5: T5-L5:

Root of Mean Square Error (RMSE): 0.185048265050244
Sum of Squared Residual: 0.719100068360212
Correlation Coef. (R): 0.999502208306672
R-Square: 0.999004664409913

Parameter	Best Estimate
----------	-------------
a1	17.7699852315831
a2	8.00653745200654
a3	0.42767787817109
a4	0.133811896038511

[Scott MacKenzie]

Am I missing something? I don't see any code in your comment.

[Alex Sha]

The below is the code for T5-L5:

Variable x,y;
Function y=A1*tanh(exp((x-A2)/A3))^A4;
Data;
x=[2:0.5:7.5 7.7 8:0.5:10.5 11.5 12];
y=[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];

when running in 1stOpt, the result will be:

Root of Mean Square Error (RMSE): 0.185048265050244
Sum of Squared Residual: 0.719100068360216
Correlation Coef. (R): 0.999502208325959
R-Square: 0.999004664448469

Parameter	Best Estimate
----------	-------------
a1	17.7700127061749
a2	8.00656197680192
a3	0.427679189766609
a4	0.133811270063284

[Image Analyst]

OK, but we don't have 1stOpt - I don't even know what it is. Is there some code that runs in MATLAB? With standard fitting functions like exist in base MATLAB or the Stats Toolbox?

[Sulaymon Eshkabilov]

It seems that 1stOpt is a 3rd party package or application. It would be nice if @Alex Sha could elaborate more on 1stOpt that demonstrates to be quite accurate. Thanks.

[Alex Sha]

More intrduction about 1stOpt could be viewd from: 7d-soft.com/en/

[Mushtaq Al-Jubbori]

My previous model is y=A1/(1+exp(A2-A3x))^1/A4, I wish to developement the function coth(x)-1/x

[Image Analyst]

@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]

Thank you, but when run the cod appear.. Undefined function or variable 'FitCoth'

Answer 2

[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.

[Image Analyst]

What is the Vandermonde model? Do you have any formula for it? If so we might be able to use the Statistics and Machine Learning Toolbox, which I believe is more common than the Curve Fitting Toolbox. Or you could pick some other of various sigmoid functions, such as the Michaelis-Menton rate equation. It would be good to know the theory behind what generated the data so we could pick a proper theoretical equation.

[Sulaymon Eshkabilov]

Vandermonde matrix can be built based on many different formulations. Its beauty and advantage is that it is computationally very efficient. E.g: Fit model: ax^2+bx+c

% y is [Measured DATA];

x = linspace(0, 5); x =x';

V = [x.^2 x ones(size(x))];

Fit_Model = V\y

[John D'Errico]

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]

@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.

Answer 3

[Image Analyst]

@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

[Mushtaq Al-Jubbori]

Very thanks for youre answer but I need fitt all curve especialy the saturation of the curve, at point the curve must be constant as in imeage

[Image Analyst]

Why??? Both Scott and I suggested a piecewise fit would be more accurate. Sure you could fit a super high order polynomial to it but it would oscillate wildly between the training points - I'm sure you've seen those demos in your classes. So give us some definite rationale for why it's absolutely necessary to have a single formula. You know that a model is not for recreating your training points -- you already know those values -- it's for predicting points that are not in your training set. So why ruin those predictions by having a single formula that won't predict them well.

[Mushtaq Al-Jubbori]

Very thanks for you and Scott, The data is relation between Track length of alpha pariclese and etchin time there fore, the diff of the function must be verify the attached imge there for, I used the function A1*tanh(U)^A4, where U=exp((x-A2)/A3), then I find the rlation between A's and energy of alpha particle.

I think polyfit doesn't work. I need one function and I find relation between parameters of this function with alpha energy

[Mushtaq Al-Jubbori]

Very thanks for you and Scott, The data is relation between Track length of alpha pariclese and etchin time there fore, the diff of the function must be verify the attached imge there for, I used the function A1*tanh(U)^A4, where U=exp((x-A2)/A3), then I find the rlation between A's and energy of alpha particle.

I think polyfit doesn't work. I need one function and I find relation between parameters of this function with alpha energy

[Image Analyst]

Well I think you could put that model into my code for fitnlm() instead of the model I used. Did you do that and try it? How well did it work out.

Otherwise, how about the model Alex used in his answer?

[Mushtaq Al-Jubbori]

I search an another equation for this I will use Alex's model

[Image Analyst]

I didn't understand that. So you're going to search for another model? Or you're going to use Alex's model? Please clarify.

Also clarify if you're going to program up your function

"the function A1*tanh(U)^A4, where U=exp((x-A2)/A3)"

Do want to use that equation (instead of another one you're going to search for, or mine, or Alex's model)? If you still want to use your tanh model, let me know if you're going to do it or you need us to do it.

Bottom line, it's not clear if you're going to search for another model, or if you're going to use one mentioned here (either Alex's, mine, Scott's, or your tanh model). What model are you're going to use?

[Mushtaq Al-Jubbori]

May be the function coth(x)-1/x good

[Image Analyst]

Is there a theory that says that your physical process should theoretically follow a coth(x)-1/x curve?

Or is it just something that looks like it fits?

[Mushtaq Al-Jubbori]

If this function ok I can explien the physical meaning

[Image Analyst]

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]

Did we can develobment the function by + or - any parameters?

Answer 4

[Scott MacKenzie]

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);

Answer 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)])

[Alex Sha]

The fitting function of "fun = p(1)*(1+exp(q(1)+q(2)*x) + p(2)*x).^-0.65;" is near same as I post out previously, that is "y = p1/(1+Exp(p2+p3*x)+p4*x)^p5" , unless put p5 as constant of "0.65"

[Image Analyst]

How did either of you come up with this model? Is it a well known model? Did you just guess? Any theory behind it? Is it a known model for his physical process of "Track length of alpha pariclese and etchin time"?

[Alex Sha]

Hi, the model function was found out automatically by 1stOpt, no any physicall meanning.

[Mushtaq Al-Jubbori]

 If the function of coth ok I expline the physical meaning of parameters

[Image Analyst]

@Mushtaq Al-Jubbori, I've seen you twice now say that you will explain the theoretical basis for using the coth(x) -1/x function, but you have not done so yet. Why not? Well never mind, I've lost interest. Good luck with your project.

I've given you code to fit that function. Only you can say it it's OK to you and meets your needs.

[Mushtaq Al-Jubbori]

Very thanks but function coth not ok, there for I cant not explin else the function is complet (no. Parameters)

[Image Analyst]

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.

Answer 6

[Mushtaq Al-Jubbori]

Thanks

Sulaymon Eshkabilov

But the cod dont run in Matlab-2016, the function optimvar undefine