浅谈matlab多变量拟合

首先申明本人是土木专业的，因为有需要要用到matlab中的拟合用途，今天好好学习了一些关于matlab多变量拟合的东西，从网上下载了一些程序，也运行了一下，就举一些实例，附上源程序吧，主要是两个自变量和三个自变量，一个因变量的拟合。让自己也更清楚，以后用起来也方便。原理就是给出一个自变量和因变量的矩阵，然后给出一个自己认为的带有未知数的拟合方程，然后付一组初始值，根据matlab返回的初始值和误差在附一组初始值，知道最后的相关系数较大，也就是误差较小时，就能拟合的比较好，写出拟合后的方程了。1.广义线性回归拟合和源码(两个自变量，一个因变量，非线性拟合）【例】这里有这样一组数据，涉及三个变量：p，t和z，要拟合出z=f(p,t)的关系式（非线性的）。zp0.811.2t609.7387520.7536.598712013.572529.632550.9387518018.9787536.5987580.13875240207512538.2212590.92530022.05544.58104.7725为了使得回归分析的结果更加直观，我调用regstats函数，编写了一个更为实用的函数：reglm，代码如下（代码中有调用方法和例子）。首先写一个M文件：functionstats=reglm(y,X,model,varnames)%多重线性回归分析或广义线性回归分析%%reglm(y,X)，产生线性回归分析的方差分析表和参数估计结果，并以表格形式显示在屏幕上.参%数X是自变量观测值矩阵，它是n行p列的矩阵.y是因变量观测值向量，它是n行1列的列向量.%%stats=reglm(y,X)，还返回一个包括了回归分析的所有诊断统计量的结构体变量stats.%%stats=reglm(y,X,model)，用可选的model参数来控制回归模型的类型.model是一个字符串，%其可用的字符串如下%'linear'带有常数项的线性模型（默认情况）%'interaction'带有常数项、线性项和交叉项的模型%'quadratic'带有常数项、线性项、交叉项和平方项的模型%'purequadratic'带有常数项、线性项和平方项的模型%%stats=reglm(y,X,model,varnames)，用可选的varnames参数指定变量标签.varnames%可以是字符矩阵或字符串元胞数组，它的每行的字符或每个元胞的字符串是一个变量的标签，它的行%数或元胞数应与X的列数相同.默认情况下，用X1,X2,…作为变量标签.%%例:%x=[215250180250180215180215250215215%136.5136.5136.5138.5139.5138.5140.5140.5140.5138.5138.5]';%y=[6.27.54.85.14.64.62.83.14.34.94.1]';%reglm(y,x,'quadratic')%%----------------------------------方差分析表----------------------------------%方差来源自由度平方和均方F值p值%回归5.000015.02773.00557.61220.0219%残差5.00001.97420.3948%总计10.000017.0018%%均方根误差(RootMSE)0.6284判定系数(R-Square)0.8839%因变量均值(DependentMean)4.7273调整的判定系数(AdjR-Sq)0.7678%%-----------------------------------参数估计-----------------------------------%变量估计值标准误t值p值%常数项30.94282011.11170.01540.9883%X10.70400.64051.09920.3218%X2-0.848729.1537-0.02910.9779%X1*X2-0.00580.0044-1.31320.2461%X1*X10.00030.00030.83840.4400%X2*X20.00520.10550.04920.9626%%Copyright2009-2010xiezhh.%$Revision:1.0.0.0$$Date:2009/12/2221:41:00$ifnargin2error('至少需要两个输入参数');endp=size(X,2);%X的列数，即变量个数ifnargin3||isempty(model)model='linear';%model参数的默认值end%生成变量标签varnamesifnargin4||isempty(varnames)varname1=strcat({'X'},num2str([1:p]'));varnames=makevarnames(varname1,model);%默认的变量标签elseifischar(varnames)varname1=cellstr(varnames);elseifiscell(varnames)varname1=varnames(:);elseerror('varnames必须是字符矩阵或字符串元胞数组');endifsize(varname1,1)~=perror('变量标签数与X的列数不一致');elsevarnames=makevarnames(varname1,model);%指定的变量标签endendST=regstats(y,X,model);%调用regstats函数进行线性回归分析，返回结构体变量STf=ST.fstat;%F检验相关结果t=ST.tstat;%t检验相关结果%显示方差分析表fprintf('\n');fprintf('------------------------------方差分析表------------------------------');fprintf('\n');fprintf('%s%7s%15s%15s%15s%12s','方差来源','自由度','平方和','均方','F值','p值');fprintf('\n');fmt='%s%13.4f%17.4f%17.4f%16.4f%12.4f';fprintf(fmt,'回归',f.dfr,f.ssr,f.ssr/f.dfr,f.f,f.pval);fprintf('\n');fmt='%s%13.4f%17.4f%17.4f';fprintf(fmt,'残差',f.dfe,f.sse,f.sse/f.dfe);fprintf('\n');fmt='%s%13.4f%17.4f';fprintf(fmt,'总计',f.dfe+f.dfr,f.sse+f.ssr);fprintf('\n');fprintf('\n');%显示判定系数等统计量fmt='%22s%15.4f%25s%10.4f';fprintf(fmt,'均方根误差(RootMSE)',sqrt(ST.mse),'判定系数(R-Square)',ST.rsquare);fprintf('\n');fprintf(fmt,'因变量均值(DependentMean)',mean(y),'调整的判定系数(AdjR-Sq)',...ST.adjrsquare);fprintf('\n');fprintf('\n');%显示参数估计及t检验相关结果fprintf('-------------------------------参数估计-------------------------------');fprintf('\n');fprintf('%8s%18s%15s%15s%12s','变量','估计值','标准误','t值','p值');fprintf('\n');fori=1:size(t.beta,1)ifi==1fmt='%8s%20.4f%17.4f%17.4f%12.4f\n';fprintf(fmt,'常数项',t.beta(i),t.se(i),t.t(i),t.pval(i));elsefmt='%10s%20.4f%17.4f%17.4f%12.4f\n';fprintf(fmt,varnames{i-1},t.beta(i),t.se(i),t.t(i),t.pval(i));endendifnargout==1stats=ST;%返回一个包括了回归分析的所有诊断统计量的结构体变量end%-----------------子函数-----------------------functionvarnames=makevarnames(varname1,model)%生成指定模型的变量标签p=size(varname1,1);varname2=[];fori=1:p-1varname2=[varname2;strcat(varname1(i),'*',varname1(i+1:end))];endvarname3=strcat(varname1,'*',varname1);switchmodelcase'linear'varnames=varname1;case'interaction'varnames=[varname1;varname2];case'quadratic'varnames=[varname1;varname2;varname3];case'purequadratic'varnames=[varname1;varname3];end调用自编的reglm函数进行二次拟合，命令如下：z=[9.7387520.7536.5987513.572529.632550.9387518.9787536.5987580.1387520.7512538.2212590.92522.05544.58104.7725];[p,t]=meshgrid([0.811.2],[60:60:300]);stats=reglm(z(:),[p(:),t(:)],'quadratic',{'p','t'});[pnew,tnew]=meshgrid(linspace(0.8,1.2,20),linspace(60,300,20));pp=pnew(:);tt=tnew(:);zhat=[ones(400,1)ppttpp.*ttpp.^2tt.^2]*stats.beta;mesh(pnew,tnew,reshape(zhat,[20,20]));holdonplot3(p(:),t(:),z(:),'k*')拟合结果：------------------------------------方差分析表------------------------------------方差来源自由度平方和均方F值p值回归5.000011548.91472309.782993.47390.0000残差9.0000222.394224.7105总计14.000011771.3089均方根误差(RootMSE)4.9710判定系数(R-Square)0.9811因变量均值(DependentMean)41.2168调整的判定系数(AdjR-Sq)0.9706-----------------------------------参数估计-----------------------------------变量估计值标准误t值p值常数项242.618869.04393.51400.0066p-513.7781137.3777-3.73990.0046t-0.36370.1212-3.00020.0150p*t0.60220.09266.50100.0001p*p272.262568.06773.99990.0031t*t-0.00030.0002-1.19460.26282.三个自变量，一个因变量clear,clcx1=[333.1