# ndm13002 / curve_fitting

contents

ndm13002 committed Apr 13, 2017
1 parent a4f6ef5 commit e15de78c005e16dfdcd26a0e85a3630a85de50ff
Showing with 131 additions and 0 deletions.
1. +50 −0 boussinesq_lookup.m
2. +21 −0 cost_logistic.m
3. +26 −0 initializeData.m
4. +8 −0 least_squares.m
5. plot1.png
6. +23 −0 problem_1_data.m
7. +3 −0 setdefaults.m
 @@ -0,0 +1,50 @@ function [sigma_z]=boussinesq_lookup(q,a,b,z) fmn= [0.1,0.02926,0.03007,0.03058 0.2,0.05733,0.05894,0.05994 0.3,0.08323,0.08561,0.08709 0.4,0.10631,0.10941,0.11135 0.5,0.12626,0.13003,0.13241 0.6,0.14309,0.14749,0.15027 0.7,0.15703,0.16199,0.16515 0.8,0.16843,0.17389,0.17739]; m = a/z; n = b/z; if (n>1.5) n=1.6; elseif (1.3<=n) if (n<=1.5) n=1.4; end elseif 1.3<=n n=1.2; end c = ones(1,4); d = ones(1,4); for i=1:4 [~,p]=min(abs(m-fmn(:,1))); M =fmn(p,1); fmn(p,1)=0; if n==1.2 t=fmn(p,2); elseif n==1.4 t=fmn(p,3); elseif n==1.6 t=fmn(p,4); end d(i) = M; c(i) = t; end b1=c(1); b2=(c(2)-c(1))/(d(2)-d(1)); b3=(((c(3)-c(2))/(d(3)-d(2)))-((c(2)-c(1))/(d(2)-d(1)))/(d(3)-d(1))); b4=(((c(4)-c(3))/(d(4)-d(3)))-((c(3)-c(2))/(d(3)-d(2)))-((c(2)-c(1))/(d(2)-d(1))))/(d(4)-d(1)); f3=b1+(b2*(m-d(1)))+(b3*(m-d(1))*(m-d(2)))+(b4*(m-d(1))*(m-d(2))*(m-d(3))); sigma_z=q*f3; end
 @@ -0,0 +1,21 @@ function [cost, gradient] = cost_logistic(a, x, y) cost = 0; gradient = 0; t = a(1)+a(2).*x; sigma = 1./(1+exp(-t)); cost = sum(-y.*log(sigma)- (1-y).*log(1-sigma)); costFun = @ (a) sum(-y.*log((1./(1+exp(-(a(1)+a(2).*x)))))-(1-y).*log(1-(1./(1+exp(-(a(1)+a(2).*x)))))); gradient = (1/length(x))*sum((sigma-y).*t); ai = [0 0]; % Set options for fminunc options = optimset('GradObj','on','MaxIter',400); % Run fminunc to obtain the optimal theta % This function will return theta and the cost [theta, cost] = fminunc(costFun, ai); t = theta(1)+theta(2).*x; sigma = 1./(1+exp(-t)); plot(x,y,'xb', x, sigma); title('Regression') xlabel('Temp (Degrees F)') ylabel('Pass or Fail (1 or 0)') end
 @@ -0,0 +1,26 @@ c = [1.0000 53.0000 1.0000; 2.0000 57.0000 1.0000; 3.0000 58.0000 1.0000; 4.0000 63.0000 1.0000; 5.0000 66.0000 0; 6.0000 66.8000 0; 7.0000 67.0000 0; 8.0000 67.2000 0; 9.0000 68.0000 0; 10.0000 69.0000 0; 11.0000 69.8000 1.0000; 12.0000 69.8000 0; 13.0000 70.2000 1.0000; 14.0000 70.2000 0; 15.0000 72.0000 0; 16.0000 73.0000 0; 17.0000 75.0000 0; 18.0000 75.0000 1.0000; 19.0000 75.8000 0; 20.0000 76.2000 0; 21.0000 78.0000 0; 22.0000 79.0000 0; 23.0000 81.0000 0]; a = c(:,1); x = c(:,2); y = c(:,3);
 @@ -0,0 +1,8 @@ function [a,fx,r2] = least_squares(Z,y) a = Z\y; Sr = sum((y-Z*a).^2); r2 = 1-Sr/sum((y-mean(y)).^2); x = Z(:,2); fx = a(1)+(a(2)*x)+(a(3)*x.^2); end
BIN +20.1 KB plot1.png
Binary file not shown.
 @@ -0,0 +1,23 @@ clear clc % part a xa=[1 2 3 4 5]'; ya=[2.2 2.8 3.6 4.5 5.5]'; Za = [ones(size(xa)), xa, xa.^-1]; [a1,fx1,r2_1] = least_squares(Za,ya); % part b xb=[3 4 5 7 8 9 11 12]'; yb=[1.6 3.6 4.4 3.4 2.2 2.8 3.8 4.6]'; Zb = [ones(size(xb)), xb, xb.^2, xb.^3]; [a2,fx2,r2_2] = least_squares(Zb,yb); % part c xc=[0.5 1 2 3 4 5 6 7 9]'; yc=[6 4.4 3.2 2.7 2.2 1.9 1.7 1.4 1.1]'; Zc = [exp(-1.5*xc),exp(-0.3*xc), exp(-0.05*xc)]; [a3,fx3,r2_3] = least_squares(Zc,yc);
 @@ -0,0 +1,3 @@ set(0, 'defaultAxesFontSize', 16) set(0,'defaultTextFontSize',14) set(0,'defaultLineLineWidth',3)