How do I perform a sin^2 or a cos^2 fit for my data?

I have the following x data and y data:

x = [0,10,20,30,40,50,60,70,80,90,100,110,120,130,140,150,160,170,180]

y = [880,1200,1900,2700,3300,3600,3000,2100,1300,860,1200,1900,2700,3150,3300,2900,2100,1300,860]

How do I perform a cos^2 (x) or a sin^2 (x) fit for the above data?

Answer 1

[Star Strider]

Perhaps this —

x = [0,10,20,30,40,50,60,70,80,90,100,110,120,130,140,150,160,170,180];
y = [880,1200,1900,2700,3300,3600,3000,2100,1300,860,1200,1900,2700,3150,3300,2900,2100,1300,860];

LvLocs = islocalmax(y, 'MinProminence',mean(y));
xlocs = x(LvLocs)

xlocs = 1x2

50 140

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

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

freq = 1/mean(diff(xlocs))

freq = 0.0111

[lb,ub] = bounds(y)

lb = 860

ub = 3600

fcn = @(b,x) b(1).*sin(2*pi*x*b(2)).^2 + b(3);

[B,fv] = fminsearch(@(b)norm(y - fcn(b,x)), rand(3,1))

B = 3x1

1.0e+03 * 2.5314 0.0006 0.9046

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

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

fv = 371.2053

txt = sprintf('$f(x) = %.2f \\cdot sin(2 \\pi x %.3f)^2 %+.3f$', B)

txt = '$f(x) = 2531.45 \cdot sin(2 \pi x 0.595)^2 +904.574$'

figure
plot(x, y)
hold on
plot(x, fcn(B,x), '-r')
hold off
grid
xlabel('X')
ylabel('Y')
text(50, 3750, txt, 'Interpreter','LaTeX')

.

[Akshay Arvind]

Thanks a lot!

[Star Strider]

As always, my pleasure!

[Akshay Arvind]

Is there a way with which I can make this fit smoother? The trend of the fit looks quite okay, but the fit is a bit crude in terms of look.

[Star Strider]

Sure!

That only requires creating a smoother ‘x’ vector with higher resolution. I also added a phase term, and tweaked the initial parameter estimates to make the fit a bit more reliable.

x = [0,10,20,30,40,50,60,70,80,90,100,110,120,130,140,150,160,170,180];
y = [880,1200,1900,2700,3300,3600,3000,2100,1300,860,1200,1900,2700,3150,3300,2900,2100,1300,860];

LvLocs = islocalmax(y, 'MinProminence',mean(y));
xlocs = x(LvLocs);
freq = 1/(2*mean(diff(xlocs)));
[lb,ub] = bounds(y);

% fcn = @(b,x) b(1).*sin(2*pi*x*b(2)).^2 + b(3);
fcn = @(b,x) b(1).*sin(2*pi*x*b(2) + b(3)).^2 + b(4);
B0 = [ub-lb; freq; rand; lb];
[B,fv] = fminsearch(@(b)norm(y - fcn(b,x)), B0)

B = 4x1

1.0e+03 * 2.5246 0.0000 -0.0000 0.9213

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

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

fv = 332.3104

mdl = fitnlm(x, y, fcn, B)                      % Calculate Statistics (Nor Required)

mdl = 
Nonlinear regression model:
    y ~ b1*sin(2*pi*x*b2 + b3)^2 + b4

Estimated Coefficients:
          Estimate         SE         tStat       pValue  
          _________    __________    _______    __________

    b1       2524.6         54.67      46.18    1.3773e-17
    b2    0.0055511    3.8132e-05     145.58    4.7733e-25
    b3    -0.046424      0.024564    -1.8899      0.078249
    b4       921.31        34.028     27.075    3.7727e-14


Number of observations: 19, Error degrees of freedom: 15
Root Mean Squared Error: 85.8
R-Squared: 0.993,  Adjusted R-Squared 0.992
F-statistic vs. constant model: 717, p-value = 2.06e-16

% txt = sprintf('$f(x) = %.2f \\cdot sin(2 \\pi x %.3f)^2 %+.3f$', B);
txt = sprintf('$f(x) = %.2f \\cdot sin(2 \\pi \\cdot %.5f \\cdot x %+.3f)^2 %+.3f$', B);

xv = linspace(min(x), max(x), 250);

figure
plot(x, y, '.:b', 'DisplayName','Data')
hold on
plot(xv, fcn(B,xv), '-r', 'DisplayName','Regression')
hold off
grid
xlabel('X')
ylabel('Y')
legend('Location','best')
title(txt, 'Interpreter','LaTeX')

The fit is better statistically as well.

.

[Akshay Arvind]

I see, thanks very much! I appreciate the help.

[Star Strider]

As always, my pleasure!