version 1.3.0.0 (47.1 KB) by
Mark Mikofski

Fit polynomial to data, forcing y-intercept to zero or arb. value and slope to zero or arb. value

[update 2014-02-19: fix error in polyfitBM, add root and select-terms tools]

This submission contains four convenience polynomial fitting functions similar to POLYFIT.

1. POLYFITZERO - fit polynomial to data, forcing y-intercept to zero.

2. POLYFITB - force y-intercept to "b".

3. POLYFITB0 - force y-intercept to "b" and slope at (0,b) to zero.

4. POLYFITBM - force y-intercept to "b" and slope at x=0 to "m", e.g.: dy/dx = m.

5. POLYFITBROOT - force intercept and root

6. POLYFITBMROOT - force intercept, slope and root

7. POLYFITBMROOTTERMS - force intercept, slope, root and terms

If you get to POLYFITBMROOTTERMS please just use POLYFITN by John D'Errico:

http://www.mathworks.com/matlabcentral/fileexchange/34765-polyfitn

Forcing the y-intercept to zero is accomplished by noting for polynomial p = [aN, ..., a3, a2, a1, a0],

IE: y = aN*y^N + ... + a3*y^3 + a2*y^2 + a1*y + a0

when x is zero, then y is the constant term, "a0". Therefore a0 = 0, or in the general case when forcing the y-intercept to an arbitrary value, a0 = b.

Forcing the slope at x=0, is accomplished by noticing that the derivative of p

IE: dy/dx = N*aN*y^(N-1) + ... + 3*a3*y^2 + 2*a2*y + a1

evaluated at x = 0 yields "a1". Therefore a1 = 0 for zero slope, or in the general case when forcing the slope to an arbitrary value, a1 = m.

Mark Mikofski (2021). polyfitZero (https://www.mathworks.com/matlabcentral/fileexchange/35401-polyfitzero), MATLAB Central File Exchange. Retrieved .

Created with
R2011b

Compatible with any release

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

Start Hunting!Create scripts with code, output, and formatted text in a single executable document.

Macarena SantillanThank you for this! Is there any way to know the gof when using this function?

Jan UlbrichCould someone explain what "root" does? I didnt find any explanation or example for that.

christos ziavrasAyca AltaySukhwinder SinghReece Nixon-LukeFantastic! Much faster than polyfit() for my data set too.

ValentinoThomas GaubertSimon StevnsDana-Adriana BotesteanuFri AniGood effort with the curve fitting function. I tried using it over an extended range for a very complex equation; its performance was below par when compared with the curve fitting function in Microsoft Excel.

J GNever mind, figured it out. Great function!

J GI am trying to use polyfitZero function to fit polynomial to data, forcing y-intercept to zero, i.e. get the polynomial to go through (2,0), however this doesn't seem to work and I'm unsure what I'm doing wrong. Thanks in advance!

a = [2.0044 2.0056 2.0021 2.0021 2.0048 2.0026 2.0035 2.0013 2.0035 2.0026];

b = [0.1006 0.0848 0.0502 0.0502 0.0909 0.0385 0.0732 0.0732 0.0896 0.0772];

scatter(a, b, 6);

hold on

p = polyfitZero(a,b,1);

f = polyval(p,a);

plot(a,f,'Color',[0.7500 0.7500 0.7500],'linewidth',1.5)

box on;

ylim([0 0.11]);

Andreas J.Lucia Morales-RivasSubhashree MishraZoegigaMark MikofskiJohn, Sorry for my slow response. You can accomplish this the same way I zeroed out the last coefficient.

E.g. to force the 1st order coefficient to zero do the following:

% ...

x = x(:);

y = y(:);

% since you want might want to include the constant term,

% change "z" to be a matrix of ones, and make it one column wider

% to include the 0th-order term

% z = zeros(dim,degree); % comment out this line which doesn't include constant

z = ones(dim,degree+1); % add this line which now includes 0th-order term

for n = 1:degree

z(:,n) = x.^(degree-n+1);

end

% you don't have to calculate the constant term,

% because it will be a column of ones (i.e. x^0 = 1)

% now remove any columns that you don't want to fit

% e.g. we want to force the 1st order coefficient to zero

orderOfZeroCoeffs = 1; % a vector of terms to omit from regression, e.g. [1,0],

% you could make "orderOfZeroCoeffs" an input to the function

% now, remove the terms whose coeffs you want to zero

orderOfNonZeroCoeffs= ~ismember(degree:-1:0, orderOfZeroCoeffs);

z = z(:,orderOfNonZeroCoeffs); % replace z with the columns removed

% there are a 1000 other ways to write this code!

% comment out the following lines

% because now it's more complicated

% pZero = z\y;

% pZero = [pZero;0]';

% instead create a vector of zeros

% any terms not calculated are already set to zero

pZero = zeros(1,degree+1);

% substitute the coefficients directly into "pZero" by index assignment

pZero(orderOfNonZeroCoeffs) = z\y;

% ...

Try that out. I'm curious what you're using this for. Also check out the other polyfit0 (fex272) function that Jennifer Sanders found above. You can modify it in exactly the same way I show here, but you will also have to modify the other output parameters.

JohnI wonder if you could modify the routine so that higher powers could also be forced to be zero. For example, the current routine set the constant coefficient to zero. But, can the 1st order coefficient also be forced to be zero. And so on...

Jennifer SandersGahh! It didn't post my initial message which was asking how to generate the structure S as you can ordinarily so with polyfit. Thus allowing you to calculate the error on the gradient. And then I found a file that does just that on the file exchange (See above post: which it posted twice?!).

Jennifer SandersOh, I just found there is an m-file that does what I want; polyfit0.m, found at:

http://www.mathworks.com/matlabcentral/fileexchange/272

Thanks anyway,

Jen

Jennifer SandersOh, I just found there is an m-file that does what I want; polyfit0.m, found at:

http://www.mathworks.com/matlabcentral/fileexchange/272

Thanks anyway,

Jen

Mark MikofskiThanks John - You are absolutely right! I originally wrote polyfitZero for a coworker who wanted to constrain the y-intercept to zero as in Excel, but later when I realized that I could use conv(), I liked it for the silly reasons that I rarely get to divide polynomials or use convolution. But you were right to point out the dangers of dividing by zero, and I should have been more explicit in the comment above that polyfitZero_v2 is not suited for functions that cross the y-axis, for all the reasons you pointed out. I wish I could delete that comment. Thanks again.

John D'ErricoSigh. I don't like to do this, but Mark's comment is simply a terribly bad idea. It might be correct in the absence of error. It might make sense if you have no points where x is actually near zero. But it is a terrible idea in any other case. Please, don't follow his advice!!!!!

What happens when you do have an exact zero in x? Yes, you get a divide by zero. But even for a near zero in any element of x, what happens is you are essentially doing a weighted regression, where the point with a near zero gets a huge weight. This will be the only point that matters in the regression analysis.

The fact is, a polynomial fit with no constant term is trivial to accomplish by a variety of means. This code uses a reasonable method, that treats the data properly and avoids risk of divide by zeros for normal data.

Again, please avoid following Mark's advice.

Mark MikofskiThere is also an even easier way to do this. Just divide by x first and then call polyfit, then use conv to multiply the 2 polynomials

function pZero2 = polyfitZero_v2(x,y,n)

q = y./x;

qnozero = q(x~=0);xnozero = x(x~=0);

pZero2 = polyfit(xnozero,qnozero,n-1);

pZero2 = conv([1 0],pZero2);

end