Polyfit not giving expected answers on rotated data

Hello, I am trying to determine the angle of movement of a target from radar data. I'm trying to use a line of best fit in simulated east-north data (with noise added). The issue came when I got different answers from the same data, but rotated 45 degrees. I wrote a script to see if polyfit was the problem:

rotang = 80;
data1 = [1 :20: 100;zeros(1,5)];
rot = [cosd(rotang), -sind(rotang); sind(rotang), cosd(rotang)];
data1_noise = zeros(2,5);
data2_noise = zeros(2,5);
for m = 1 : 5
    data1_noise(1,m) = data1(1,m) + randn*5;
    data1_noise(2,m) = data1(2,m) + randn*5;
    data2_noise(:,m) = rot * data1_noise(:,m) ;
end
data2 = [(1 :20: 100) *cosd(rotang); (1 :20: 100) *sind(rotang)];
p1 = polyfit(data1_noise(1,:),data1_noise(2,:),1);
p2 = polyfit(data2_noise(1,:),data2_noise(2,:),1);
y1 = polyval(p1,data1_noise(1,:));
y2 = polyval(p2,data2_noise(1,:));
azimuth1 = atand(p1(1));
azimuth2 = atand(p2(1));
disp(azimuth1);disp(azimuth2-rotang);
figure(1);clf;hold on
plot(data1(1,:),data1(2,:),'b')
plot(data2(1,:),data2(2,:),'r')
plot(data1_noise(1,:),data1_noise(2,:),'xb')
plot(data2_noise(1,:),data2_noise(2,:),'xr')
plot(data1_noise(1,:),y1,'--b')
plot(data2_noise(1,:),y2,'--r')
axis equal
hold off

I would expect azimuth2 to be azimuth1 + the rotation angle, since data2 is just data1 rotated, but it is not. Depending on the random noise, it can be completely different:

Does anyone understand what I am doing wrong?

Many thanks

[Torsten]

5*randn is quite a lot of noise ... Orthogonal regression instead of usual least-squares regression should give you the expected relation between p1(1) and p2(1) - independent of the noise introduced.

https://en.wikipedia.org/wiki/Deming_regression

[Joe Bailey]

Deming worked, thanks

Answer 1

[Cris LaPierre]

polyfit finds the best fit in a least-squares sense. That means it uses the vertical distance between each point and the fit line to determine the best fit. The vertical distance changes as you rotate the line. See the visualizaions below from the curve fitter app. Note the change in scale of the x axis between the two; data1_noise is on the top, and data2_noise is on the bottom.

[Cris LaPierre]

What you want to do is instead use perpendicular distance. Then you would get what you expect - the fit line angle only differs by rotang.

There may be a simpler way to do this, but I followed the example given here: Fitting an Orthogonal Regression using PCA

% create your data
rotang = 80;
data1 = [1 :20: 100;zeros(1,5)];
rot = [cosd(rotang), -sind(rotang); sind(rotang), cosd(rotang)]

rot = 2×2

0.1736 -0.9848 0.9848 0.1736

<mw-icon class=""></mw-icon>

<mw-icon class=""></mw-icon>

data1_noise = data1 + randn(size(data1))*5;
data2_noise = rot * data1_noise;

% Compute az of data1_noise
az1 = computeAz(data1_noise')

az1 = 0.2322

% Compute az of data2_noise
az2 = computeAz(data2_noise')

●

az2 = 80.2322

az2-rotang

ans = 0.2322

function az = computeAz(X)
[n,p] = size(X);
meanX = mean(X,1);
[coeff,score,roots] = pca(X);
dirVect = coeff(:,1);

Xfit1 = repmat(meanX,n,1) + score(:,1)*coeff(:,1)';

t = [min(score(:,1))-.2, max(score(:,1))+.2];
endpts = [meanX + t(1)*dirVect'; meanX + t(2)*dirVect'];
plot(endpts(:,1),endpts(:,2),'k-');

X1 = [X(:,1) Xfit1(:,1)];
X2 = [X(:,2) Xfit1(:,2)];

hold on
plot(X1',X2','b-', X(:,1),X(:,2),'bo');
hold off
axis equal
grid on

delta = diff(endpts);
slope = delta(2)/delta(1);
az = atand(slope);
end

[Joe Bailey]

Thanks