No BSD License  

Highlights from
Regression Analysis

from Regression Analysis by Gerhard Nieuwenhuiys
Do regression analysis (curve fitting) on 2D, 3D graphs and Differential equations.

regres(fstring,v,xexp,yexp);
% Author : Gerhard Nieuwenhuiys
% Date of last revision : 2009/02/03
% Thank you for downloading / using regres V1. Should you have any inquiries
% about operating instruction or you would like to share your comments for
% improvement, feel free to contact me at gerhardn654@yahoo.com
% The function regres(fstring,v,xexp,yexp) requires the following 4 inputs.
% regres can only handle explicit functions of the form y=f(x).
% fstring = The model to be fitted to the data points.
% v = The values where fminsearch should start searching for a solution.
% The number of elements in v should be the same as the number of independent
% parameters in fstring. v is only able to accept a max of 10 searchable
% parameters.
% xexp = experimental values corresponding to the independent variable
% yexp = experimental values corresponding to the dependent variable
%
% Example 1.
% regres('p1.*x+p2',[1 1],[0 1 2 3 4],[4 5 6 7 8]);
%  p1=1 p2=4 R2=1

function regres(fstring,v,xexp,yexp);
[w1,w2,w3,w4,w5,w6,w7,w8,w9,w10]=deal(strfind(fstring,'p1'),strfind(fstring,'p2'),strfind(fstring,'p3'),...
    strfind(fstring,'p4'),strfind(fstring,'p5'),strfind(fstring,'p6'),strfind(fstring,'p7'),strfind(fstring,'p8'),...
    strfind(fstring,'p9'),strfind(fstring,'p10'));
[k1,k2,k3,k4,k5,k6,k7,k8,k9,k10]=deal(length(w1),length(w2),length(w3),length(w4),length(w5),length(w6),length(w7),...
    length(w8),length(w9),length(w10));
[x,fval,exitflag,output]=fminsearch(@errfun,v,[],fstring,xexp,yexp);
sse=fval;
ys=sum(yexp)/length(yexp+xexp);
sst=sum((yexp-ys).^2);
R2=1-(sse/sst);
x1=min(xexp):(max(xexp)-min(xexp))/1000:max(xexp);
if k1>=1 & k2>=1 & k3>=1 & k4>=1 & k5>=1 & k6>=1 & k7>=1 & k8>=1 & k9>=1 & k10>=1;
    disp([' p1=',num2str(x(1)),' p2=',num2str(x(2)),' p3=',num2str(x(3)),' p4=',num2str(x(4)),' p5=',num2str(x(5)),...
        ' p6=',num2str(x(6)),' p7=',num2str(x(7)),' p8=',num2str(x(8)),' p9=',num2str(x(9)),' p10=',num2str(x(10)),...
        ' R2=',num2str(R2)]);
    f=inline(fstring,'x','p1','p2','p3','p4','p5','p6','p7','p8','p9','p10');
    y1=f(x1,x(1),x(2),x(3),x(4),x(5),x(6),x(7),x(8),x(9),x(10));
elseif k1>=1 & k2>=1 & k3>=1 & k4>=1 & k5>=1 & k6>=1 & k7>=1 & k8>=1 & k9>=1;
    disp([' p1=',num2str(x(1)),' p2=',num2str(x(2)),' p3=',num2str(x(3)),' p4=',num2str(x(4)),' p5=',num2str(x(5)),...
        ' p6=',num2str(x(6)),' p7=',num2str(x(7)),' p8=',num2str(x(8)),' p9=',num2str(x(9)),' R2=',num2str(R2)]);
    f=inline(fstring,'x','p1','p2','p3','p4','p5','p6','p7','p8','p9');
    y1=f(x1,x(1),x(2),x(3),x(4),x(5),x(6),x(7),x(8),x(9));
elseif k1>=1 & k2>=1 & k3>=1 & k4>=1 & k5>=1 & k6>=1 & k7>=1 & k8>=1;
    disp([' p1=',num2str(x(1)),' p2=',num2str(x(2)),' p3=',num2str(x(3)),' p4=',num2str(x(4)),' p5=',num2str(x(5)),...
        ' p6=',num2str(x(6)),' p7=',num2str(x(7)),' p8=',num2str(x(8)),' R2=',num2str(R2)]);
    f=inline(fstring,'x','p1','p2','p3','p4','p5','p6','p7','p8');
    y1=f(x1,x(1),x(2),x(3),x(4),x(5),x(6),x(7),x(8));
elseif k1>=1 & k2>=1 & k3>=1 & k4>=1 & k5>=1 & k6>=1 & k7>=1;
    disp([' p1=',num2str(x(1)),' p2=',num2str(x(2)),' p3=',num2str(x(3)),' p4=',num2str(x(4)),' p5=',num2str(x(5)),...
        ' p6=',num2str(x(6)),' p7=',num2str(x(7)),' R2=',num2str(R2)]);
    f=inline(fstring,'x','p1','p2','p3','p4','p5','p6','p7');
    y1=f(x1,x(1),x(2),x(3),x(4),x(5),x(6),x(7));
elseif k1>=1 & k2>=1 & k3>=1 & k4>=1 & k5>=1 & k6>=1;
    disp([' p1=',num2str(x(1)),' p2=',num2str(x(2)),' p3=',num2str(x(3)),' p4=',num2str(x(4)),' p5=',num2str(x(5)),...
        ' p6=',num2str(x(6)),' R2=',num2str(R2)]);
    f=inline(fstring,'x','p1','p2','p3','p4','p5','p6');
    y1=f(x1,x(1),x(2),x(3),x(4),x(5),x(6));
elseif k1>=1 & k2>=1 & k3>=1 & k4>=1 & k5>=1;
    disp([' p1=',num2str(x(1)),' p2=',num2str(x(2)),' p3=',num2str(x(3)),' p4=',num2str(x(4)),' p5=',num2str(x(5)),...
        ' R2=',num2str(R2)]);
    f=inline(fstring,'x','p1','p2','p3','p4','p5');
    y1=f(x1,x(1),x(2),x(3),x(4),x(5));
elseif k1>=1 & k2>=1 & k3>=1 & k4>=1;
    disp([' p1=',num2str(x(1)),' p2=',num2str(x(2)),' p3=',num2str(x(3)),' p4=',num2str(x(4)),' R2=',num2str(R2)]);
    f=inline(fstring,'x','p1','p2','p3','p4');
    y1=f(x1,x(1),x(2),x(3),x(4));
elseif k1>=1 & k2>=1 & k3>=1;
    disp([' p1=',num2str(x(1)),' p2=',num2str(x(2)),' p3=',num2str(x(3)),' R2=',num2str(R2)]);
    f=inline(fstring,'x','p1','p2','p3');
    y1=f(x1,x(1),x(2),x(3));
elseif k1>=1 & k2>=1;
    disp([' p1=',num2str(x(1)),' p2=',num2str(x(2)),' R2=',num2str(R2)]);
    f=inline(fstring,'x','p1','p2');
    y1=f(x1,x(1),x(2));
elseif k1>=1;
    disp([' p1=',num2str(x(1)),' R2=',num2str(R2)]);
    f=inline(fstring,'x','p1');
    y1=f(x1,x(1));
end;
plot(xexp,yexp,'r.',x1,y1,'b-'),axis image,grid on;

function err=errfun(x,fstring,xexp,yexp);
[w1,w2,w3,w4,w5,w6,w7,w8,w9,w10]=deal(strfind(fstring,'p1'),strfind(fstring,'p2'),strfind(fstring,'p3'),...
    strfind(fstring,'p4'),strfind(fstring,'p5'),strfind(fstring,'p6'),strfind(fstring,'p7'),strfind(fstring,'p8'),...
    strfind(fstring,'p9'),strfind(fstring,'p10'));
[k1,k2,k3,k4,k5,k6,k7,k8,k9,k10]=deal(length(w1),length(w2),length(w3),length(w4),length(w5),length(w6),length(w7),...
    length(w8),length(w9),length(w10));
if k1>=1 & k2>=1 & k3>=1 & k4>=1 & k5>=1 & k6>=1 & k7>=1 & k8>=1 & k9>=1 & k10>=1;
    [p1,p2,p3,p4,p5,p6,p7,p8,p9,p10]=deal(x(1),x(2),x(3),x(4),x(5),x(6),x(7),x(8),x(9),x(10));
    f=inline(fstring,'x','p1','p2','p3','p4','p5','p6','p7','p8','p9','p10');
    g=f(xexp,p1,p2,p3,p4,p5,p6,p7,p8,p9,p10);
elseif k1>=1 & k2>=1 & k3>=1 & k4>=1 & k5>=1 & k6>=1 & k7>=1 & k8>=1 & k9>=1;
    [p1,p2,p3,p4,p5,p6,p7,p8,p9]=deal(x(1),x(2),x(3),x(4),x(5),x(6),x(7),x(8),x(9));
    f=inline(fstring,'x','p1','p2','p3','p4','p5','p6','p7','p8','p9');
    g=f(xexp,p1,p2,p3,p4,p5,p6,p7,p8,p9);
elseif k1>=1 & k2>=1 & k3>=1 & k4>=1 & k5>=1 & k6>=1 & k7>=1 & k8>=1;
    [p1,p2,p3,p4,p5,p6,p7,p8]=deal(x(1),x(2),x(3),x(4),x(5),x(6),x(7),x(8));
    f=inline(fstring,'x','p1','p2','p3','p4','p5','p6','p7','p8');
    g=f(xexp,p1,p2,p3,p4,p5,p6,p7,p8);
elseif k1>=1 & k2>=1 & k3>=1 & k4>=1 & k5>=1 & k6>=1 & k7>=1;
    [p1,p2,p3,p4,p5,p6,p7]=deal(x(1),x(2),x(3),x(4),x(5),x(6),x(7));
    f=inline(fstring,'x','p1','p2','p3','p4','p5','p6','p7');
    g=f(xexp,p1,p2,p3,p4,p5,p6,p7);
elseif k1>=1 & k2>=1 & k3>=1 & k4>=1 & k5>=1 & k6>=1;
    [p1,p2,p3,p4,p5,p6]=deal(x(1),x(2),x(3),x(4),x(5),x(6));
    f=inline(fstring,'x','p1','p2','p3','p4','p5','p6');
    g=f(xexp,p1,p2,p3,p4,p5,p6);
elseif k1>=1 & k2>=1 & k3>=1 & k4>=1 & k5>=1;
    [p1,p2,p3,p4,p5]=deal(x(1),x(2),x(3),x(4),x(5));
    f=inline(fstring,'x','p1','p2','p3','p4','p5');
    g=f(xexp,p1,p2,p3,p4,p5);
elseif k1>=1 & k2>=1 & k3>=1 & k4>=1;
    [p1,p2,p3,p4]=deal(x(1),x(2),x(3),x(4));
    f=inline(fstring,'x','p1','p2','p3','p4');
    g=f(xexp,p1,p2,p3,p4);
elseif k1>=1 & k2>=1 & k3>=1;
    [p1,p2,p3]=deal(x(1),x(2),x(3));
    f=inline(fstring,'x','p1','p2','p3');
    g=f(xexp,p1,p2,p3);
elseif k1>=1 & k2>=1;
    [p1,p2]=deal(x(1),x(2));
    f=inline(fstring,'x','p1','p2');
    g=f(xexp,p1,p2);
elseif k1>=1;
    p1=x(1);
    f=inline(fstring,'x','p1');
    g=f(xexp,p1);
end;
y=double(g);
err=sum((y-yexp).^2);

Contact us at files@mathworks.com