Polyfitn

I wanted to compare the speed of my SURFIT with that of POLYFITN - POLYVALN. My problem is that POLYFITN apparently accepts only a vector as dependent variable, whereas I have to get a polynomial approximation of a matrix. How should I do that?

An example from my work: imagine that I have a data matrix Z (512 by 511) which I want to approximate by a polynomial of 2-nd degree in vertical direction and 0-th degree in horizontal direction. This means I want to use a bipolynomial of degrees [2 0]. With SURFIT the problem is solved easily:

B=surfit(Z,[2 0]); %B = my bipolynomial approximation of Z (regression)

What should I do with POLYFITN-POLYVALN to get B?

Yes, the comparison was done correctly. I have noticed that the ratio of errors depends on data. Occasionally, I used real data for comparison, which allowed to notice an inacceptably big difference in accuracy. Below I have included a model of a piece of these data as default.
 

[B,B2]=john;
Elapsed time is 0.063000 seconds. %POLY...N
Elapsed time is 0.000000 seconds. %SURFIT

JV = 1.1931e+005 95907 %Errors

 function [B,B2]=john(Z)
%======================
if nargin<1
y=zeros(100,160);
y(1:30,20:60)=-8;
y(80:100,40:120)=-8;
Z=y+randn(size(y));
end

[R,C]=size(Z);
[X,Y]=meshgrid((1:C)/C,(1:R)/R);
z=[X(:),Y(:)];
ZZ=Z(:);
tic,
mod=polyfitn(z,ZZ,2);
B=polyvaln(mod,z);
toc
B=reshape(B,R,C);

tic,B2=surfit(Z,2);toc
JV=sum([ZZ-B(:),ZZ-B2(:)].^2)

does what i need it to.

note that Vassili Pastushenko wrote SURFIT.

Hi all,

Good job John as usually you do. I planed to use bootstrap approach for some educational experimental works because uncertainties can touch both x and y variables. In that goal i compared your program with my own bootstrap program using polyfit. That gave similar results except a) as expected bootstrap was slower, but that is not very important in my case, b) bootstrap failed when data points were under around 10 with a polyfit warning message stating an ill-conditioned matrix. This is due to replicated points in the bootstrap samples.

Maybe i could manage bootstrap in a more correct way. Anyway your program succeed more quickly and with reasonable uncertainties even when there was noise in the independant variable x.

So, thanks a lot John. I know that if i have some clever task to do, i before have to check if John has done it.

This is excellent work. Thanks John.
I have a question. I have a series points (x, y,z), but no any equation yet. These data was obtained from test. I am not sure if I can get the function z=f(x,y) with polyfitn? Could you help me? Thank you very much.

It's very easy to use and works excellent. Has been very useful to my research. Thank you.

Well, I need to try it out first, so i gave it four star. If it works on fitting my DIC displacement distribution, I gave you the whole star you want. But after all, thank you any way..

This error is incomplete, since it does not show the calling syntax. As such, I cannot say why the problem arose. Please send me a direct e-mail with any questions, and I will resolve them immediately.

Hello, the software is excellent, but there is a very minor error, it does not check for that the modelterm matrix passed in is not empty, resulting in an mtimes error in this case. I discovered this because of my poor coding, so I suppose I shouldn't expect polyfitn to check this for me!

Many thanks! just what I needed. saved me a lot of coding

John hello,

Even running this example of yours,

n = 1000;
x = rand(n,2);
y = exp(sum(x,2)) + randn(n,1)/100;
p = polyfitn(x,y,3)

...I'm getting this error:
??? Undefined function or method 'polyfitn' for input arguments of type 'double'.

Any idea of what the problem might be???

Thanks,

Agustin

John, could you add to the polyfitn package description file mentioned in your ReadMe file or any relevant reference on multidimensional polynomial regression.

"In additions to these m-files, there are several other files of interest. For those who have an interest in understanding how polyfitn works (quite simple really) you should look to the file understanding_polyfitn.rtf, in the doc subdirectory."

The underlying theory will be very helpful for me. Thanks in advance ...

Hallo John,
Nice work!
I am a newbie in matlab and still got some questions about 3d fit.
I have a 3d dataset F(128,128,16) and would like to get the polynomial function F(x,y,z), how should I do it?

thanks!

Great work !
Thank you very much. It is really easy to use, allows a lot of user specifications and is quite fast imho.

James,

You don't need every possible combination of (x,y,z) to be known to estimate a model, certainly not a first order model like this. To estimate the 3 unknown parameters in this model, in theory, you need at least 3 data points at a minimum. In practice, you won't get very good estimates of those parameters with only 3 data points.

In the past, I often had clients ask me how much data they needed to obtain. After all, data is often scarce, expensive, hard to obtain. It costs money, and time. You can't have too much data, at least if it is good data.

In practice, if you end up with 10 or 20 points to estimate these three parameters, that will (probably) be entirely adequate. More data will give you higher quality estimates. Even 6 points, carefully chosen, if the noise is low enough, will often be adequate. Only you know the quality of your data, and how much you have, as well as how accurately you need the model to be estimated.

If you have further questions, I can only suggest you send me an e-mail, with your data, and any other information you have on the problem.

John

If you are trying to place downloaded files into the Matlab/toolbox directories, this is generally wrong. Those directories are cached by matlab when it starts up, so any files that you put in those directories will not be seen until you restart matlab. Also, anytime you make any changes to files that you put in those directories, the changes will not be recognized by matlab.

Better is to create your own new (and separate) directory for downloaded files. Add that directory on your search path, using the pathtool function, or using addpath and then savepath.

As I explained the first time, you need to put polyfitn into a directory that is on your search path.

help addpath
help savepath
help pathtool

Your particular test problem gives me the results:

uv = rand(100,2);
w = sin(sum(uv,2));
p = polyfitn(uv,w,'u, v');

polyn2sympoly(p)
ans =
    0.76152*u + 0.73572*v

Nothing special needs be done except putting it in a directory where matlab can find it.

Toby - you found a bug in the model parser. Repaired now. Thanks, John

Great function, really easy to use! I am not sure how it compares to the new surface fitting toolbox from Matlab (R2009).

If you have data that falls on a straight line in three dimensions, then you need to use a tool designed to solve that problem. Polyfitn is not that tool.

I presume that you wish to find the line that passes through these points, assuming that each variable has noise on it. This is called a total least squares problem. The solution is commonly formulated using a singular value decomposition. That is, assume that your data is in an nx3 array, X, where each row of the array X corresponds to a single data point in three dimensions.

A line in N dimensions (here N = 3) is defined by the equation

   P(t) = P0 + P1*t

where P0 and P1 are vectors of length 3. Any point on that line is given by some value of the parameter t. This would be called a parametric form of the line. It has the advantage that is any "slopes" of the line might be zero or infinite, we never have a singularity.

How can we compute the vectors P0 and P1? P0 is simple. We simply compute the mean of our data.

   P0 = mean(X,1);

Note that P0 is not unique, as ANY point on the line would suffice to define a line. To estimate the line equation using the svd, the mean of your data is best. P1 is also easy enough to generate. It is the right singular vector that corresponds to the maximum singular value of the matrix with P0 subtracted off. Lets try it out...

Here is some data that falls along a straight line in three dimensions.

n = 100;
X = bsxfun(@plus,rand(n,1)*[2 3 5],randn(1,3)) + randn(n,3)/25;
plot3(X(:,1),X(:,2),X(:,3),'.')

P0 = mean(X,1)
P0 =
       2.3547 0.11682 0.85517

[U,S,V] = svd(bsxfun(@minus,X,P0),0);

diag(S)
ans =
       18.253
      0.40999
      0.35826

The straight line is a good fit when the maximum singular value is significantly larger than the other singular values. Here, that maximum singular value is roughly 44 times as large as the next larges singular value. The right singular vectors are the columns of V.

V
V =
        0.322 0.87967 -0.34999
        0.489 -0.47108 -0.73414
      0.81068 -0.065251 0.58185

The corresponding right singular vector is the first column of V. I'll transpose it to be consistent in shape with P0.

P1 = V(:,1).'
P1 =
        0.322 0.489 0.81068

Is P1 close to being correct? It turns out that we can scale P1 by any constant and it will be equally valid. I'll scale it here so the first element is 2. If we look back at the way I generated the data, it would originally have been [2 3 5] had there been no noise applied. How well did I do?

P1*(2/P1(1))
ans =
            2 3.0373 5.0353

I'd say that was pretty close.