Curve Fitting with erfc

27 次查看（过去 30 天）
Russell Evans
Russell Evans 2020-3-26
评论： Walter Roberson 2020-3-27
Greetings,
I am trying to fit the attached data (Ydata vs time) with this function attached as a picture.
The coefficients for fitting are B and D and must be positive. L is 1000, Q is -0.0038
I am trying to use cftool with the custom equation
((D*-0.0038/1000^2)*B) .*exp((B.^2) .*(((D.*x))./1000.^2)) .*erfc(B.*(sqrt(((D.*x))./1000^2)))
but it just says my imputs are not real. It would be helpful to see how to fit the data.
  2 个评论
John D'Errico
John D'Errico 2020-3-26
编辑：John D'Errico 2020-3-26
Why does your custom equation not look like the model you claim to want to use?
Here, for example, we see
exp((B.^2) .*((D.*x)./1000.^2))
what appears to be a term like D*t/L^2.
But then in the next term, if x is supposed to be t, then all you have is sqrt(t).
As far as knowing why your inputs are not real, we cannot conjecture what other wrong things you did, since you already have make a last one significant mistake. That suggests your skills at coding are not incredibly great. While there is nothing overtly wrong with that, it does suggest that you have surely made at least a few other mistakes in your code, especially given the vague error you have described. If you want help, the easiest way to get help is to show what you tried. Then someone can just correct your mistakes.
As well, when you have an error, rather than just trying to describe what you think is pertinent about the error message, show the complete error text, so everything in red.
Russell Evans
Russell Evans 2020-3-26
I edited the question to reflect
((D*-0.0038/1000^2)*B) .*exp((B.^2) .*(((D.*x))./1000.^2)) .*erfc(B.*(sqrt(((D.*x))./1000^2)))
This inital typo was my fault. "x" is supposed to be "t"
As you mention, my skills at coding are not incredibly great. Thanks for pointing that out.
No red error shows. I do get a yellow warning message
Error in fittype expression ==> ((D.*-0.0038./1000.^2).*B) .*exp((B.^2) .*(((D.*x))./1000.^2)) .*erfc(B.*sqrt(x))
??? Input must be real and full.
I have made sure that my x data is real even when imput into erfc(1/sqrt(x)) outside of cftool

请先登录，再进行评论。

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

采纳的回答

Walter Roberson
Walter Roberson 2020-3-26
I was able to track this down.
What is happening is that cftool likes to draw the smooth line, from before the start of the x data, to after the end of the x data. It uses a 5% margin, starting the projection at min(x) - (max(x)-min(x))*0.05 . That happens to be negative for your data, and the sqrt(x) goes complex and erfc complains.
Unfortunately there are no controls for where the line is drawn. When I changed the sqrt(x) to sqrt(max(x,0)) to prevent square root of negatives, the smooth line doesn't appear and the fit is not good.
The easiest solution, if you have the symbolic toolbox, is build the function for sum-of-squared residues, and minimize the function:
syms D B
x = time; y = ypap;
residue = sum((((D*-0.0038/1000^2)*B) .*exp((B.^2) .*((D.*x)./1000.^2)) .*erfc(B.*sqrt(x)) - y).^2);
F = matlabFunction(residue, 'vars', {[B D]});
goodBD = fminsearch(F, [1 2]);
goodB = goodBD(1);
goodD = goodBD(2);
The fminsearch is fast, and the solution it gets you is pretty much within round-off noise of the best solution I have been able to find with more careful analysis.
  4 个评论
显示 2更早的评论
Alex Sha
Alex Sha 2020-3-27
Hi, refer the result below (with the function: y=((D*(-0.0038)/1000^2)*B)*exp((B^2)*(((D*x))/1000^2))*erfc(B*(sqrt(((D*x))/1000^2)))), the result looks good but be careful in the end of chart lines (with the blue circle), it seems to have a problem with that fitting function.
Root of Mean Square Error (RMSE): 1.09533665147734E-8
Sum of Squared Residual: 4.07919209223658E-15
Correlation Coef. (R): 0.99976230572428
R-Square: 0.999524667947129
Parameter Best Estimate
---------- -------------
b 122.725546738343
d 5.47654549122705
Walter Roberson
Walter Roberson 2020-3-27
Why did you set minimum for B and D to be 1, when earlier you just said "positive" ?
syms D B real
x = time; y = ypap;
ERFC = @(G)simplify(arrayfun(@(H)piecewise(imag(H)~=0,1e10,erfc(H)),G));
I = @(X) ((D*-0.0038/1000^2)*B) .*exp((B.^2) .*(((D.*X))./1000.^2)) .*ERFC(B.*(sqrt(((D.*X))./1000^2)));
parts = arrayfun(I, x)-y;
residue = simplify(sum(simplify(parts.^2)));
F = matlabFunction(residue, 'vars', {[B D]}, 'file', 'residue_fun.m');
goodBD = fminsearch(F, [100 1]);
goodB = goodBD(1);
goodD = goodBD(2);
The starting point was important for fast convergence.

请先登录，再进行评论。

更多回答（0 个）

类别

在 Help Center 和 File Exchange 中查找有关 Get Started with Curve Fitting Toolbox 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by