Clear Filters
Clear Filters

How do I reorder a regression plot in the desired sequence?

3 views (last 30 days)
I am working on a set of datapoints like so:
x = [270 280 290 300 310 320 330 340 350 0 10 20 30 40 50 60 70 80 90]
y = [10000 9000 5500 2500 900 2500 5500 9000 10000 9000 5500 2500 900 2500 5500 9000 10000]
I have defined a sine^2(x) function with a phase shift of 45 degree, for fitting along these data points, because my data points are shifted that way.
I define x1 = 1:numel(x);
Using set(gca,'xTick',x1,'XTickLabel',x), I know that I can display in the x-axis order 270,...,0,...,90, without which I would get them in the x-axis order 0...90,...,270,...,350.
But, how do I apply that order for the fitting function itself? That doesn't seem to work.
I attach the code for your reference.
function[]=plotdata(filename)
S = load(filename);
C = struct2cell(S);
M = cell2mat(C);
x = M(:,1)
y = M(:,end)
x1 = 1:numel(x);
mean(y)
xlocs = [270 0 90]
%freq = 1/(2*mean(diff(xlocs))) %diff(xlocs)
freq = 1 / (2*mean(diff(xlocs)))
[lb,ub] = bounds(y)
fcn = @(b,x)b(1).*cos(2*pi*x*b(2)+b(3)+(pi/4)).^2+b(4)
B0 = [ub-lb; freq; 0; lb]
myfun = @(b)norm(fcn(b,x) - y);
[B,fv] = fminsearch(myfun,B0)
xv = linspace(max(x),min(x),1000); %A smoother x vector
figure
plot(x, y, '*', 'DisplayName','Data')
hold on
plot(xv, fcn(B,xv), '-r', 'DisplayName','Regression')
hold off
%set(gca,'xTick',x1,'XTickLabel',x)
end
  2 Comments
Mathieu NOE
Mathieu NOE on 2 Jul 2024
why is it important for you to have the data sorted this way ?
that has no influence on the fitting process
Akshay Arvind
Akshay Arvind on 2 Jul 2024
Hello Mathieu, while I also realized in principle the fitting is absolutely fine without the sorting of the data, it would be something nice if I can have it sorted (just more presentability of the plot). It also reflects that I use a polarizer from 270 degree to 90 degree.

Sign in to comment.

Accepted Answer

Mathieu NOE
Mathieu NOE on 2 Jul 2024
hello again
you could do that
NB that the cos^2 is not what I used
A simple sin (or cos if you prefer works better - your y data does not have a squared cos shape
also I added a few missing y points in your data (x and y vectors must be same length)
to make the math a bit simpler for me , I prefer to have x monotonic, so for the 0 to 90 deg segment I added 360 deg
and I simply retrieved again this 360 deg shift to x1 to have the xticklabels as you wanted
hope it helps !
% x = [270 280 290 300 310 320 330 340 350 0 10 20 30 40 50 60 70 80 90];
x = [[270 280 290 300 310 320 330 340 350] [0 10 20 30 40 50 60 70 80 90]+360];
% y = [10000 9000 5500 2500 900 2500 5500 9000 10000 9000 5500 2500 900 2500 5500 9000 10000];
y = [10000 9000 5500 2500 900 2500 5500 9000 10000 9000 5500 2500 900 2500 5500 9000 10000 9000 5500 ];
%%%%%%%%%%%%% main code %%%%%%%%%%%%%%%%%
ym = mean(y); % Estimate offset
yu = max(y);
yl = min(y);
yr = (yu-yl)/2; % Estimate Amplitude
% Find zero-crossings
zt = find_zc(x,y,ym);
per = mean(diff(zt)); % Estimate period
fre = 1/per; % Estimate FREQUENCY
phe = -2*pi*zt(1)/per; % Estimate phase
% stationnary sinus fit
fit = @(b,t) b(1)*(sin(2*pi*t*b(2) + b(3))) + b(4); % Objective Function to fit
fcn = @(b) norm(fit(b,x) - y); % Least-Squares cost function
sol = fminsearch(fcn, [yr; fre; phe; ym;]); % Minimise Least-Squares
sol(3) = mod(sol(3),2*pi);
yp = fit(sol,x);
% display sinus fit parameters in command window
amplitude = sol(1);
frequency_Hz = sol(2)
frequency_Hz = 0.0125
phase_offset_rad = sol(3)
phase_offset_rad = 5.4708
DC_offset = sol(4)
DC_offset = 5.6044e+03
% plot
figure(1),
plot(x, y, 'b',zt, ym*ones(size(zt)), 'dr',x, yp, '-g')
legend('data','zc points','model fit');
xlim([x(1) x(end)]);
x1 = x;
ind = x1>=360;
x1(ind) = x1(ind) - 360;
set(gca,'xTick',x,'XTickLabel',num2str(x1'))
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [Zx] = find_zc(x,y,threshold)
% positive slope "zero" crossing detection, using linear interpolation
y = y - threshold;
zci = @(data) find(diff(sign(data))>0); %define function: returns indices of +ZCs
ix=zci(y); %find indices of + zero crossings of x
ZeroX = @(x0,y0,x1,y1) x0 - (y0.*(x0 - x1))./(y0 - y1); % Interpolated x value for Zero-Crossing
Zx = ZeroX(x(ix),y(ix),x(ix+1),y(ix+1));
end

More Answers (0)

Products


Release

R2022b

Community Treasure Hunt

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

Start Hunting!