浅谈matlab多变量拟合
matlab回归（拟合）总结（一元、多元）
matlab 回归（拟合）总结前言1、学三条命令polyfit(x,y,n)---拟合成一元幂函数（一元多次） regress(y,x)----可以多元，nlinfit(x,y,’fun ’,beta0) (可用于任何类型的函数，任意多元函数，应用范围最主，最万能的)2、同一个问题，这三条命令都可以使用，但结果肯定是不同的，因为拟合的近似结果，没有唯一的标准的答案。
相当于咨询多个专家。
3、回归的操作步骤：根据图形（实际点），选配一条恰当的函数形式（类型）---需要数学理论与基础和经验。
（并写出该函数表达式的一般形式，含待定系数）------选用某条回归命令求出所有的待定系数。
所以可以说，回归就是求待定系数的过程（需确定函数的形式）一、回归命令一元多次拟合polyfit(x,y,n);一元回归polyfit;多元回归regress---nlinfit(非线性)二、多元回归分析对于多元线性回归模型(其实可以是非线性，它通用性极高)： e x x y pp ++++=βββ 110设变量12,,,p x x x y 的n 组观测值为12(,,,)1,2,,i i ip i x x x y i n =记 ⎪⎪⎪⎪⎪⎭⎫⎝⎛=np n n p p x x x x x x x x x x 212222111211111，⎪⎪⎪⎪⎪⎭⎫ ⎝⎛=n y y y y 21，则⎪⎪⎪⎪⎪⎭⎫⎝⎛=p ββββ 10 的估计值为排列方式与线性代数中的线性方程组相同（），拟合成多元函数---regress使用格式：左边用b=[b, bint, r, rint, stats]右边用=regress(y, x)或regress(y, x, alpha) ---命令中是先y 后x,---须构造好矩阵x(x 中的每列与目标函数的一项对应) ---并且x 要在最前面额外添加全1列/对应于常数项---y 必须是列向量---结果是从常数项开始---与polyfit 的不同。
matlab拟合工具箱拟合方法
matlab拟合工具箱拟合方法Matlab拟合工具箱是Matlab软件中的一个功能强大的工具箱，它提供了多种拟合方法，用于拟合数据集并找到最佳的拟合曲线。
本文将介绍Matlab拟合工具箱的几种常用的拟合方法。
一、线性拟合（Linear Fit）线性拟合是最简单和最常用的拟合方法之一。
线性拟合假设拟合曲线为一条直线，通过最小二乘法求解最佳拟合直线的斜率和截距。
线性拟合可以用于解决一些简单的线性关系问题，例如求解两个变量之间的线性关系、求解直线运动的速度等。
二、多项式拟合（Polynomial Fit）多项式拟合是一种常见的拟合方法，它假设拟合曲线为一个多项式函数。
多项式拟合可以适用于一些非线性的数据集，通过增加多项式的阶数，可以更好地拟合数据。
在Matlab拟合工具箱中，可以通过设置多项式的阶数来进行多项式拟合。
三、指数拟合（Exponential Fit）指数拟合是一种常用的非线性拟合方法，它假设拟合曲线为一个指数函数。
指数拟合可以用于拟合一些呈指数增长或指数衰减的数据集。
在Matlab拟合工具箱中，可以使用指数拟合函数来拟合数据集，并得到最佳的拟合曲线参数。
四、对数拟合（Logarithmic Fit）对数拟合是一种常见的非线性拟合方法，它假设拟合曲线为一个对数函数。
对数拟合可以用于拟合一些呈对数增长或对数衰减的数据集。
在Matlab拟合工具箱中，可以使用对数拟合函数来拟合数据集，并得到最佳的拟合曲线参数。
五、幂函数拟合（Power Fit）幂函数拟合是一种常用的非线性拟合方法，它假设拟合曲线为一个幂函数。
幂函数拟合可以用于拟合一些呈幂函数增长或幂函数衰减的数据集。
在Matlab拟合工具箱中，可以使用幂函数拟合函数来拟合数据集，并得到最佳的拟合曲线参数。
六、指数幂函数拟合（Exponential Power Fit）指数幂函数拟合是一种常见的非线性拟合方法，它假设拟合曲线为一个指数幂函数。
指数幂函数拟合可以用于拟合一些呈指数幂函数增长或指数幂函数衰减的数据集。
matlab中多项式拟合方法
MATLAB中多项式拟合方法一、概述在科学计算和工程领域，多项式拟合是一种常用的数据拟合方法。
MATLAB作为一种强大的数学计算软件，提供了多种多项式拟合的函数和工具，可以方便地进行数据拟合和分析。
二、多项式拟合的原理多项式拟合是利用多项式函数来拟合已知的数据点，使得多项式函数与实际数据点的残差最小化。
多项式函数可以表达为:\[ y(x) = a_0 + a_1x + a_2x^2 + ... + a_nx^n \]其中，\(y(x)\)为拟合函数，\(a_0, a_1, a_2,...,a_n\)为多项式系数，\(x\)为自变量。
拟合的目标是通过确定系数的取值，使得多项式函数和实际数据点的误差最小。
三、MATLAB中的多项式拟合函数MATLAB提供了多种函数和工具来进行多项式拟合，常用的函数包括polyfit、polyval和polyfitn等。
1. polyfit函数polyfit函数用于多项式拟合，其调用格式为:\[ p = polyfit(x, y, n) \]其中，\(x\)为自变量数据，\(y\)为因变量数据，\(n\)为拟合的多项式阶数。
函数返回一个多项式系数向量\(p\)，可以使用polyval函数计算拟合的多项式函数值。
2. polyval函数polyval函数用于计算多项式函数的值，其调用格式为:\[ y_fit = polyval(p, x) \]其中，\(p\)为多项式系数向量，\(x\)为自变量数据，\(y_fit\)为拟合的多项式函数值。
3. polyfitn函数polyfitn函数是MATLAB中的一个拟合工具箱，可以进行更复杂的多项式拟合和数据分析，包括多变量多项式拟合、非线性多项式拟合等。
四、多项式拟合的应用多项式拟合在科学研究和工程实践中有着广泛的应用，例如数据分析、曲线拟合、信号处理等领域。
1. 数据分析多项式拟合可用于分析实验数据，拟合实验结果，从而得出数据之间的关系和规律。
Matlab多变量二次多项式拟合
一、对Y 总做线性多项式拟合：0112288......Y b b X b X b X =+++ 设置显著性水平为0.05，拟合得到：B=[ 0b ,1b ,………., 8b ]= [-60.0349 12.5809 2.2002 -12.9863 20.4145 0.0266 5.1430 17.2416 151.6779]对应的置信区间为： -161.4058 41.3359 -7.5870 32.7488 -25.5706 29.9709 -33.5089 7.5362 -0.3096 41.1386 -2.5989 2.6520 0.9830 9.3030 -3.2810 37.7642 -64.0209 367.3767r 2= 0.7454 (越接近于1，回归效果越显著)，F= 2.5616， p= 0.1163，（p>0.05, 可知回归模型不成立）。
残差图如下：从残差图可以看出，除第一个数据和最后一个数据的残差离零点均较远，说明这两个数据可视为异常点，去掉这两个数据之后再做拟合得到：B=[ 0b ,1b ,………., 8b ]= [-478.8 15.7 1.8 -85.343 2.8 24.7 135.3 1131.9]对应的置信区间为：-1048.791.1 7.5 23.9 -8 11.6 -183.5 12.8 10.5 75.5 -1.1 6.7 -2 51.4 -25.8 296.4 -206.7 2470.4r 2= 0.9690 (越接近于1，回归效果越显著)，F= 19.5530， p= 0.0023，（p<0.05, 可知回归模型成立）。
残差图如下：从残差图可以看出，数据的残差离零点均较近，且残差的置信区间均包含零点，这说明回归模型能较好的符合原始数据。
预测值与实测值的比较： Y 总 预测值 Y 总 实测值 相对误差42.7458 46.48 -8.03% 40.5008 41.99 -3.55% 44.9358 43.18 4.07% 66.1358 66.06 0.11% 42.9108 42.3 1.44% 35.76 35.76 0.00% 64.4508 61.67 4.51% 40.9608 38.18 7.28% 52.46 52.46 0.00% 62.5008 61.89 0.99% 39.2758 39.2 0.19% 60.4758 58.72 2.99% 64.9108 66.4 -2.24% 62.665866.4-5.62%从上表可以看出，预测值和实测值的误差都在10%以内，说明该拟合模型能很好的预测实验值。
浅谈matlab多变量拟合
首先申明本人是土木专业的，因为有需要要用到matlab中的拟合用途，今天好好学习了一些关于matlab多变量拟合的东西，从网上下载了一些程序，也运行了一下，就举一些实例，附上源程序吧，主要是两个自变量和三个自变量，一个因变量的拟合。
让自己也更清楚，以后用起来也方便。
原理就是给出一个自变量和因变量的矩阵，然后给出一个自己认为的带有未知数的拟合方程，然后付一组初始值，根据matlab返回的初始值和误差在附一组初始值，知道最后的相关系数较大，也就是误差较小时，就能拟合的比拟好，写出拟合后的方程了。
1.广义线性回归拟合和源码(两个自变量，一个因变量，非线性拟合〕【例】这里有这样一组数据，涉及三个变量：p，t 和z，要拟合出z = f(p,t) 的关系式〔非线性的〕。
z p 0.8 1 1.2t60 9.73875 20.75 36.5987120 13.5725 29.6325 50.93875180 18.97875 36.59875 80.13875240 2075125 38.22125 90.925300 22.055 44.58 104.7725为了使得回归分析的结果更加直观，我调用regstats函数，编写了一个更为实用的函数：reglm，代码如下〔代码中有调用方法和例子〕。
首先写一个M文件：function stats = 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 = [215 250 180 250 180 215 180 215 250 215215% 136.5 136.5 136.5 138.5 139.5 138.5 140.5 140.5 140.5 138.5 138.5]';% y = [6.2 7.5 4.8 5.1 4.6 4.6 2.8 3.1 4.3 4.9 4.1]';% reglm(y,x,'quadratic')%% ----------------------------------方差分析表----------------------------------% 方差来源自由度平方和均方F值p值% 回归 5.0000 15.0277 3.0055 7.6122 0.0219% 残差 5.0000 1.9742 0.3948% 总计10.0000 17.0018%% 均方根误差(Root MSE) 0.6284 判定系数(R-Square) 0.8839% 因变量均值(Dependent Mean) 4.7273 调整的判定系数(Adj R-Sq) 0.7678%% -----------------------------------参数估计-----------------------------------% 变量估计值标准误t值p 值% 常数项30.9428 2011.1117 0.0154 0.9883% X1 0.7040 0.640 5 1.0992 0.3218% X2 -0.8487 29.1537 -0.0291 0.9779%X1*X2 -0.0058 0.0044 -1.3132 0.2461% X1*X1 0.00030.0003 0.8384 0.4400% X2*X2 0.0052 0.1055 0.0492 0.9626 %% Copyright 2009 - 2010 xiezhh. % $Revision: 1.0.0.0 $ $Date: 2009/12/22 21:41:00 $if nargin < 2 error('至少需要两个输入参数');endp = size(X,2); % X的列数，即变量个数if nargin < 3 || isempty(model) model = 'linear'; % model参数的默认值end% 生成变量标签varnamesif nargin < 4 || isempty(varnames) varname1 =strcat({'X'},num2str([1:p]')); varnames = makevarnames(varname1,model); % 默认的变量标签else if ischar(varnames) varname1 = cellstr(varnames); elseifiscell(varnames) varname1 = varnames(:); else error('varnames 必须是字符矩阵或字符串元胞数组'); end if size(varname1,1) ~= p error('变量标签数与X的列数不一致'); else varnames = 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,'均方根误差(Root MSE)',sqrt(ST.mse),'判定系数(R-Square)',ST.rsquare);fprintf('\n');fprintf(fmt,'因变量均值(Dependent Mean)',mean(y),'调整的判定系数(AdjR-Sq)',... ST.adjrsquare); fprintf('\n');fprintf('\n');% 显示参数估计及t检验相关结果fprintf('-------------------------------参数估计-------------------------------');fprintf('\n');fprintf('%8s%18s%15s%15s%12s','变量','估计值','标准误','t 值','p值');fprintf('\n');for i = 1:size(t.beta,1) if i == 1 fmt ='%8s%20.4f%17.4f%17.4f%12.4f\n'; fprintf(fmt,'常数项',t.beta(i),t.se(i),t.t(i),t.pval(i)); else fmt ='%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)); endendif nargout == 1 stats = ST; % 返回一个包括了回归分析的所有诊断统计量的构造体变量end% -----------------子函数-----------------------function varnames = makevarnames(varname1,model)% 生成指定模型的变量标签p = size(varname1,1);varname2 = [];for i = 1:p-1 varname2 = [varname2;strcat(varname1(i),'*',varname1(i+1:end))];endvarname3 =strcat(varname1,'*',varname1);switch model case 'linear' varnames = varname1; case'interaction' varnames = [varname1;varname2]; case 'quadratic' varnames = [varname1;varname2;varname3]; case 'purequadratic' varnames = [varname1;varname3];end 调用自编的reglm函数进展二次拟合，命令如下：>> z = [9.73875 20.75 36.5987513.5725 29.6325 50.9387518.97875 36.59875 80.1387520.75125 38.22125 90.92522.055 44.58 104.7725];>> [p,t] = meshgrid([0.8 1 1.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) pp tt pp.*tt pp.^2 tt.^2]*stats.beta;>> mesh(pnew,tnew,reshape(zhat,[20,20]));>> hold on>> plot3(p(:),t(:),z(:),'k*')拟合结果：------------------------------------方差分析表------------------------------------方差来源自由度平方和均方F值p值回归 5.0000 11548.9147 2309.7829 93.4739 0.0000残差9.0000 222.3942 24.7105总计14.0000 11771.3089 均方根误差(Root MSE) 4.9710 判定系数(R-Square) 0.9811因变量均值(DependentMean) 41.2168 调整的判定系数(Adj R-Sq) 0.9706-----------------------------------参数估计----------------------------------- 变量估计值标准误t值p值常数项242.6188 69.0439 3.5140 0.0066 p -513.7781 137.3777 -3.7399 0.0046 t -0.3637 0.1212 -3.0002 0.0150 p*t0.6022 0.0926 6.5010 0.0001 p*p 272.2625 68.0677 3.9999 0.0031 t*t -0.0003 0.0002 -1.1946 0.26282.三个自变量，一个因变量clear,clcx1=[333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 328.15 330.65 333.15 335.65 338.15 340.65 343.15 333.15 333.15 333.15 323.15 325.65 345.65 348.15]';x2=[1.19 1.206 1.228 1.23 1.252 1.27 1.277 1.31 1.35 1.39 1.43 1.23 1.23 1.23 1.23 1.23 1.2 1.2 1.2 1.2 1.2 1.26 1.26 1.26 1.26 1.26 1.23 1.23 1.23 1.23 1.23 1.23 1.23 1.15 1.47 1.51 1.23 1.23 1.23 1.23]';x3=[80 80 80 80 80 80 80 80 80 80 80 77 78 79 80 81 67 68 69 70 71 86 87 88 89 90 80 80 80 80 80 80 80 80 80 80 80 80 80 80]'; y=[59.49 55.16 50.18 49.78 45.75 42.96 41.96 37.87 33.96 30.83 28.29 47.92 48.54 49.19 49.78 50.42 47.49 48.21 48.9 49.63 50.32 47.8 48.38 48.91 49.47 50.04 50.49 50.14 49.79 49.45 49.13 48.81 48.5 74.13 26.18 24.39 51.22 50.85 48.2147.92]';X=[x1,x2,x3];ymin=min(y);y=y-ymin;fx=(b,x1,x2,x3)(b(1)+b(2)*x1+b(3)*x2+b(4)*x3+b(5)*x 1.^2+b(6)*x2.^2+b(7)*x3.^2+b(8)*x1.*x2+b(9)*x1.*x3+b(10)*x2.*x3+b(11)*x1.^3+b(12)*x2.^3+b( 13)*x3.^3)./(1+b(14)*exp(b(15)*x1+b(16)*x2+b(17)*x3+b(18)*x1.*x2+b(19)*x1.*x3+b(20)*x2.*x3)) ;fx2=(b,X,y)(b(1)+b(2)*X(:,1)+b(3)*X(:,2)+b(4)*X(:,3)+b(5)*X(:,1).^2+b(6)*X(:,2).^2+b(7)*X(:,3).^2 +b(8)*X(:,1).*X(:,2)+b(9)*X(:,1).*X(:,3)+b(10)*X(:,2).*X(:,3)+b(11)*X(:,1).^3+b(12)*X(:,2).^3+b(13) *X(:,3).^3)./(1+b(14)*exp(b(15)*X(:,1)+b(16)*X(:,2)+b(17)*X(:,3)+b(18)*X(:,1).*X(:,2)+b(19)*X(:,1). *X(:,3)+b(20)*X(:,2).*X(:,3)));bm=[105091.513651451,1328.,-711027.452435498,-1213.,-1.625,934239 .742471165,-25.43,-1301.,10.67,-642.1,0.,-244987.606559315,0.9581,9.986e-05,-0.19651,13.74,0.32436, -0.,0.7,-0.50883];b=bm;forl=1:10 b=lsqcurvefit(fx2,b,X,y); b=nlinfit(X,y,fx2,b);endbm1=mean(x1);m2=mean(x2);m3=me an(x3);r1=range(x1);r2=range(x2);r3=range(x3);ry=range(y);x1a=min(x1);x1b=max(x1);x2a=min(x2);x2b=max(x2);x3a= min(x3);x3b=max(x3);ya=min(y);yb=max(y);n=length(y);str=num2str([1:n]');figure(1),clfplot3(x1,x2,y ,'o')stem3(x1,x2,y,'filled')text(x1,x2,y+.04*ry,str,'fontsize',12)pause(.0001)holdon[x11,x22]=meshgrid(x1a:r1/75:x1b,x2a:r2/75:x2b);y1=fx(bm,x11,x22,m3);surf(x11,x22,y1)axis([x1a x1b x2a x2b ya yb])alpha(.85)shading interpaxis tightpause(1.0001)%clf% forl=1:10% plot3(x1,x2,y,'o')% stem3(x1,x2,y,'filled')% text(x1,x2,y+.04*ry,str,'fontsize',12)% pause(.0001)% holdon% m3=x3a+l*r3/10;% y1=fx(bm,x11,x22,m3);% surf(x11,x22,y1)% axis([x1a x1b x2a x2b ya yb])% alpha(.4)% shading interp% axis tight% pause(.5001)%endxlabel('X1'),ylabel('X2'),zlabel('Y')figure(2),clf[x11,x33]=meshgrid(x1a:r1/75:x1b,x3a:r3/75:x3b);plot3(x1,x3,y,'o')stem3(x1,x3,y,'filled')text(x1,x3,y+.04*ry,str,'fontsize',12)pause(.0001)holdony2=fx(bm,x11,m2,x33);surf(x11,x33,y2)axis([x1a x1b x3a x3b ya yb])alpha(.85)shading interpaxis tightpause(5.0001)%clf% forl=1:10% plot3(x1,x3,y,'o')% stem3(x1,x3,y,'filled')% text(x1,x3,y+.04*ry,str,'fontsize',12)% pause(.0001)% holdon% m2=x2a+(l-1)*r2/10;% y2=fx(bm,x11,m2,x33);% surf(x11,x33,y2)% axis([x1a x1b x3a x3b ya yb])% alpha(.4)% shading interp% axis tight% pause(.5001)%endxlabel('X1'),ylabel('X3'),zlabel('Y')figure(3),clfplot3(x2,x3,y,'o')stem3(x2,x3,y,'filled')text(x2,x3,y+.0 4*ry,str,'fontsize',12)pause(.0001)holdon[x22,x33]=meshgrid(x2a:r2/75:x2b,x3a:r3/75:x3b);y3=fx(bm,m1,x22,x33);surf(x22,x33,y3)axis([x2a x2b x3a x3b ya yb])alpha(.85)shading interpaxis tightpause(5.0001)%clf% forl=1:10% plot3(x2,x3,y,'o')% stem3(x2,x3,y,'filled')% text(x2,x3,y+.04*ry,str,'fontsize',12)% pause(.0001)% holdon% m1=x1a+(l-1)*r1/10;% y3=fx(bm,m1,x22,x33);% surf(x22,x33,y3)% axis([x2a x2b x3a x3b ya yb])% alpha(.4)% shading interp% axis tight% pause(.5001)%endxlabel('X2'),ylabel('X3'),zlabel('Y')disp([' x1 x2 x3 y yhat'])yhat=fx(b,x1,x2,x3);[x1,x2,x3,y+ymin,yhat+ymin]SSy=var(y)*(n-1)RSS=(y-yhat)'*(y-yhat)Ra qaure=(SSy-RSS)/SSy。
matlab拟合多元函数
matlab拟合多元函数
matlab拟合多元函数
Matlab可以通过"fitlm" 命令来拟合多元函数。
具体操作步骤如下：
1. 将需要拟合的函数按照一定的格式写出来，比如：
y = β1*x1 + β2*x2 + β3*x3 + ε
其中，y表示因变量，x1、x2、x3表示自变量，β1、β2、β3表示对应自变量的系数，ε表示误差。
2. 构建数据矩阵并输入到fitlm中进行拟合，比如：
data = [x1 x2 x3 y];
model = fitlm(data, 'y ~ x1 + x2 + x3');
其中，data是一个包含自变量和因变量的矩阵，model是拟合出来的模型。
3. 获取拟合结果，比如：
coefficients = model.Coefficients.Estimate;
R_squared = model.Rsquared.Ordinary;
其中，coefficients是拟合出来的系数向量，R_squared是判定系数。
总之，通过fitlm命令能够比较方便地实现多元函数的拟合。
浅论matlab多变量拟合
首先申明本人是土木专业的，因为有需要要用到matlab中的拟合用途，今天好好学习了一些关于matlab多变量拟合的东西，从网上下载了一些程序，也运行了一下，就举一些实例，附上源程序吧，主要是两个自变量和三个自变量，一个因变量的拟合。
让自己也更清楚，以后用起来也方便。
原理就是给出一个自变量和因变量的矩阵，然后给出一个自己认为的带有未知数的拟合方程，然后付一组初始值，根据matlab返回的初始值和误差在附一组初始值，知道最后的相关系数较大，也就是误差较小时，就能拟合的比较好，写出拟合后的方程了。
1.广义线性回归拟合和源码(两个自变量，一个因变量，非线性拟合）【例】这里有这样一组数据，涉及三个变量：p，t 和z，要拟合出 z = f(p,t) 的关系式（非线性的）。
z p 0.8 1 1.2t60 9.73875 20.75 36.5987120 13.5725 29.6325 50.93875180 18.97875 36.59875 80.13875240 2075125 38.22125 90.925300 22.055 44.58 104.7725为了使得回归分析的结果更加直观，我调用regstats函数，编写了一个更为实用的函数：reglm，代码如下（代码中有调用方法和例子）。
首先写一个M文件：function stats = 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 = [215 250 180 250 180 215 180 215 250 215 215% 136.5 136.5 136.5 138.5 139.5 138.5 140.5 140.5 140.5 138.5 138.5]'; % y = [6.2 7.5 4.8 5.1 4.6 4.6 2.8 3.1 4.3 4.9 4.1]';% reglm(y,x,'quadratic')%% ----------------------------------方差分析表----------------------------------% 方差来源自由度平方和均方 F值 p值% 回归 5.0000 15.0277 3.0055 7.6122 0.0219% 残差 5.0000 1.9742 0.3948% 总计 10.0000 17.0018%% 均方根误差(Root MSE) 0.6284 判定系数(R-Square) 0.8839% 因变量均值(Dependent Mean) 4.7273 调整的判定系数(AdjR-Sq) 0.7678%% -----------------------------------参数估计-----------------------------------% 变量估计值标准误 t值 p 值% 常数项 30.9428 2011.1117 0.0154 0.9883 % X1 0.7040 0.6405 1.0992 0.3218 % X2 -0.8487 29.1537 -0.0291 0.9779 % X1*X2 -0.0058 0.0044 -1.3132 0.2461 % X1*X1 0.0003 0.0003 0.8384 0.4400 % X2*X2 0.0052 0.1055 0.0492 0.9626 %% Copyright 2009 - 2010 xiezhh.% $Revision: 1.0.0.0 $ $Date: 2009/12/22 21:41:00 $if nargin < 2error('至少需要两个输入参数');endp = size(X,2); % X的列数，即变量个数if nargin < 3 || isempty(model)model = 'linear'; % model参数的默认值end% 生成变量标签varnamesif nargin < 4 || isempty(varnames)varname1 = strcat({'X'},num2str([1:p]'));varnames = makevarnames(varname1,model); % 默认的变量标签elseif ischar(varnames)varname1 = cellstr(varnames);elseif iscell(varnames)varname1 = varnames(:);elseerror('varnames 必须是字符矩阵或字符串元胞数组');endif size(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,'均方根误差(Root MSE)',sqrt(ST.mse),'判定系数(R-Square)',ST.rsquare);fprintf('\n');fprintf(fmt,'因变量均值(Dependent Mean)',mean(y),'调整的判定系数(Adj R-Sq)',...ST.adjrsquare);fprintf('\n');fprintf('\n');% 显示参数估计及t检验相关结果fprintf('-------------------------------参数估计-------------------------------');fprintf('\n');fprintf('%8s%18s%15s%15s%12s','变量','估计值','标准误','t值','p值'); fprintf('\n');for i = 1:size(t.beta,1)if i == 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)); endendif nargout == 1stats = ST; % 返回一个包括了回归分析的所有诊断统计量的结构体变量end% -----------------子函数-----------------------function varnames = makevarnames(varname1,model)% 生成指定模型的变量标签p = size(varname1,1);varname2 = [];for i = 1:p-1varname2 = [varname2;strcat(varname1(i),'*',varname1(i+1:end))]; endvarname3 = strcat(varname1,'*',varname1);switch modelcase 'linear'varnames = varname1;case 'interaction'varnames = [varname1;varname2];case 'quadratic'varnames = [varname1;varname2;varname3];case 'purequadratic'varnames = [varname1;varname3];end调用自编的reglm函数进行二次拟合，命令如下：>> z = [9.73875 20.75 36.5987513.5725 29.6325 50.9387518.97875 36.59875 80.1387520.75125 38.22125 90.92522.055 44.58 104.7725];>> [p,t] = meshgrid([0.8 1 1.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) pp tt pp.*tt pp.^2 tt.^2]*stats.beta;>> mesh(pnew,tnew,reshape(zhat,[20,20]));>> hold on>> plot3(p(:),t(:),z(:),'k*')拟合结果：------------------------------------方差分析表------------------------------------方差来源自由度平方和均方 F值 p值回归 5.0000 11548.9147 2309.7829 93.4739 0.0000 残差 9.0000 222.3942 24.7105总计 14.0000 11771.3089均方根误差(Root MSE) 4.9710 判定系数(R-Square) 0.9811因变量均值(Dependent Mean) 41.2168 调整的判定系数(AdjR-Sq) 0.9706-----------------------------------参数估计----------------------------------- 变量估计值标准误 t值 p值常数项 242.6188 69.0439 3.5140 0.0066 p -513.7781 137.3777 -3.7399 0.0046t -0.3637 0.1212 -3.0002 0.0150p*t 0.6022 0.0926 6.5010 0.0001p*p 272.2625 68.0677 3.9999 0.0031t*t -0.0003 0.0002 -1.1946 0.26282.三个自变量，一个因变量clear,clcx1=[333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 328.15 330.65 333.15 335.65 338.15 340.65 343.15 333.15 333.15 333.15 323.15 325.65 345.65 348.15]';x2=[1.19 1.206 1.228 1.23 1.252 1.27 1.277 1.31 1.35 1.39 1.43 1.23 1.23 1.23 1.23 1.23 1.2 1.2 1.2 1.2 1.2 1.26 1.26 1.26 1.26 1.26 1.23 1.23 1.23 1.23 1.23 1.23 1.23 1.15 1.47 1.51 1.23 1.23 1.23 1.23]';x3=[80 80 80 80 80 80 80 80 80 80 80 77 78 79 80 81 67 68 69 70 71 86 87 88 89 90 80 80 80 80 80 80 80 80 80 80 80 80 80 80]';y=[59.49 55.16 50.18 49.78 45.75 42.96 41.96 37.87 33.96 30.83 28.29 47.92 48.5449.19 49.78 50.42 47.49 48.21 48.9 49.63 50.32 47.8 48.38 48.91 49.47 50.04 50.4950.14 49.79 49.45 49.13 48.81 48.5 74.13 26.18 24.39 51.22 50.85 48.21 47.92]'; X=[x1,x2,x3];ymin=min(y);y=y-ymin;fx=@(b,x1,x2,x3)(b(1)+b(2)*x1+b(3)*x2+b(4)*x3+b(5)*x1.^2+b(6)*x2.^2+b(7)*x3.^2 +b(8)*x1.*x2+b(9)*x1.*x3+b(10)*x2.*x3+b(11)*x1.^3+b(12)*x2.^3+b(13)*x3.^3)./(1 +b(14)*exp(b(15)*x1+b(16)*x2+b(17)*x3+b(18)*x1.*x2+b(19)*x1.*x3+b(20)*x2.*x3)); fx2=@(b,X,y)(b(1)+b(2)*X(:,1)+b(3)*X(:,2)+b(4)*X(:,3)+b(5)*X(:,1).^2+b(6)*X(:, 2).^2+b(7)*X(:,3).^2+b(8)*X(:,1).*X(:,2)+b(9)*X(:,1).*X(:,3)+b(10)*X(:,2).*X(: ,3)+b(11)*X(:,1).^3+b(12)*X(:,2).^3+b(13)*X(:,3).^3)./(1+b(14)*exp(b(15)*X(:,1 )+b(16)*X(:,2)+b(17)*X(:,3)+b(18)*X(:,1).*X(:,2)+b(19)*X(:,1).*X(:,3)+b(20)*X( :,2).*X(:,3)));bm=[105091.513651451,1328.10332025611,-711027.452435498,-1213.61405762992,-1.8 8264106646625,934239.742471165,-25.5844409887743,-1301.90766627356,10.51891749 78167,-642.229950374061,0.00221335659769481,-244987.606559315,0.15540437371958 1,9.28886223888986e-05,-0.0142397533119651,13.4903417277274,0.0213803812532436 ,-0.00141251443766222,0.000377042917999337,-0.0845412180650883];b=bm;for l=1:10b=lsqcurvefit(fx2,b,X,y);b=nlinfit(X,y,fx2,b);endbm1=mean(x1);m2=mean(x2);m3=mean(x3);r1=range(x1); r2=range(x2);r3=range(x3);ry=range(y);x1a=min(x1);x1b=max(x1);x2a=min(x2);x2b=max(x2);x3a=min(x3);x3b=max(x3);ya=min(y);yb=max(y);n=length(y);str=num2str([1:n]');figure(1),clfplot3(x1,x2,y,'o')stem3(x1,x2,y,'filled')text(x1,x2,y+.04*ry,str,'fontsize',12)pause(.0001)hold on[x11,x22]=meshgrid(x1a:r1/75:x1b,x2a:r2/75:x2b);y1=fx(bm,x11,x22,m3);surf(x11,x22,y1)axis([x1a x1b x2a x2b ya yb])alpha(.85)shading interpaxis tightpause(1.0001)%clf% for l=1:10% plot3(x1,x2,y,'o')% stem3(x1,x2,y,'filled')% text(x1,x2,y+.04*ry,str,'fontsize',12)% pause(.0001)% hold on% m3=x3a+l*r3/10;% y1=fx(bm,x11,x22,m3);% surf(x11,x22,y1)% axis([x1a x1b x2a x2b ya yb])% alpha(.4)% shading interp% axis tight% pause(.5001)% endxlabel('X1'),ylabel('X2'),zlabel('Y')figure(2),clf[x11,x33]=meshgrid(x1a:r1/75:x1b,x3a:r3/75:x3b); plot3(x1,x3,y,'o')stem3(x1,x3,y,'filled')text(x1,x3,y+.04*ry,str,'fontsize',12)pause(.0001)hold ony2=fx(bm,x11,m2,x33);surf(x11,x33,y2)axis([x1a x1b x3a x3b ya yb])alpha(.85)shading interpaxis tightpause(5.0001)%clf% for l=1:10% plot3(x1,x3,y,'o')% stem3(x1,x3,y,'filled')% text(x1,x3,y+.04*ry,str,'fontsize',12)% pause(.0001)% hold on% m2=x2a+(l-1)*r2/10;% y2=fx(bm,x11,m2,x33);% surf(x11,x33,y2)% axis([x1a x1b x3a x3b ya yb])% alpha(.4)% shading interp% axis tight% pause(.5001)% endxlabel('X1'),ylabel('X3'),zlabel('Y')figure(3),clfplot3(x2,x3,y,'o')stem3(x2,x3,y,'filled')text(x2,x3,y+.04*ry,str,'fontsize',12)pause(.0001)hold on[x22,x33]=meshgrid(x2a:r2/75:x2b,x3a:r3/75:x3b);y3=fx(bm,m1,x22,x33);surf(x22,x33,y3)axis([x2a x2b x3a x3b ya yb])alpha(.85)shading interpaxis tightpause(5.0001)%clf% for l=1:10% plot3(x2,x3,y,'o')% stem3(x2,x3,y,'filled')% text(x2,x3,y+.04*ry,str,'fontsize',12)% pause(.0001)% hold on% m1=x1a+(l-1)*r1/10;% y3=fx(bm,m1,x22,x33);% surf(x22,x33,y3)% axis([x2a x2b x3a x3b ya yb])% alpha(.4)% shading interp% axis tight% pause(.5001)% endxlabel('X2'),ylabel('X3'),zlabel('Y')disp([' x1 x2 x3 y yh at'])yhat=fx(b,x1,x2,x3);[x1,x2,x3,y+ymin,yhat+ymin]SSy=var(y)*(n-1)RSS=(y-yhat)'*(y-yhat) Raqaure=(SSy-RSS)/SSy。
Matlab中的数据拟合方法探究
Matlab中的数据拟合方法探究引言：数据拟合是数据分析中的重要任务之一，它通过数学模型来描述和预测观测数据的变化规律。
Matlab作为一种功能强大的科学计算软件，提供了丰富的工具和函数来进行数据拟合。
本文将探究Matlab中常用的数据拟合方法，包括最小二乘法拟合、多项式拟合、曲线拟合以及局部加权回归拟合等。
一、最小二乘法拟合最小二乘法拟合是一种常用的线性拟合方法，它通过最小化观测值与拟合曲线之间的误差平方和来确定模型参数。
在Matlab中，可以使用linefit函数来进行最小二乘法拟合。
该函数可以根据观测数据点的坐标，拟合出一条直线，并返回拟合直线的斜率和截距。
二、多项式拟合多项式拟合是一种非线性拟合方法，它通过多项式函数来逼近观测数据。
在Matlab中，可以使用polyfit函数进行多项式拟合。
该函数可以拟合出一个多项式模型，并返回各个系数的值。
用户可以根据实际需求选择多项式的次数，以达到最佳拟合效果。
三、曲线拟合曲线拟合是一种更加灵活的非线性拟合方法，它通过拟合多条曲线来逼近观测数据。
在Matlab中，可以使用curvefit函数进行曲线拟合。
该函数可以根据用户提供的初始参数值，拟合出一个曲线模型，并返回各个参数的最优估计值。
曲线拟合可以适用于各种非线性的数据。
四、局部加权回归拟合局部加权回归拟合是一种非参数的非线性拟合方法，它通过引入权重函数，对不同数据点赋予不同的权重，来逼近观测数据的变化趋势。
在Matlab中，可以使用loess函数进行局部加权回归拟合。
该函数可以根据用户提供的带宽参数和权重函数，拟合出一条光滑的曲线，并返回各个数据点的拟合值。
五、综合应用与讨论在实际数据分析中，我们往往需要综合应用不同的拟合方法，以获得更加准确的结果。
例如，我们可以先使用最小二乘法拟合得到一个初步的线性模型，然后再通过多项式拟合或曲线拟合来进一步修正模型。
在这个过程中，我们还可以使用交叉验证等方法来评估模型的拟合效果。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
首先申明本人是土木专业的，因为有需要要用到matlab中的拟合用途，今天好好学习了一些关于matlab多变量拟合的东西，从网上下载了一些程序，也运行了一下，就举一些实例，附上源程序吧，主要是两个自变量和三个自变量，一个因变量的拟合。
让自己也更清楚，以后用起来也方便。
原理就是给出一个自变量和因变量的矩阵，然后给出一个自己认为的带有未知数的拟合方程，然后付一组初始值，根据matlab返回的初始值和误差在附一组初始值，知道最后的相关系数较大，也就是误差较小时，就能拟合的比较好，写出拟合后的方程了。
1.广义线性回归拟合和源码(两个自变量，一个因变量，非线性拟合）【例】这里有这样一组数据，涉及三个变量：p，t 和z，要拟合出z = f(p,t) 的关系式（非线性的）。
z p 0.8 1 1.2t60 9.73875 20.75 36.5987120 13.5725 29.6325 50.93875180 18.97875 36.59875 80.13875240 2075125 38.22125 90.925300 22.055 44.58 104.7725为了使得回归分析的结果更加直观，我调用regstats函数，编写了一个更为实用的函数：reglm，代码如下（代码中有调用方法和例子）。
首先写一个M文件：function stats = 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 = [215 250 180 250 180 215 180 215 250 215 215% 136.5 136.5 136.5 138.5 139.5 138.5 140.5 140.5 140.5 138.5 138.5]'; % y = [6.2 7.5 4.8 5.1 4.6 4.6 2.8 3.1 4.3 4.9 4.1]';% reglm(y,x,'quadratic')%% ----------------------------------方差分析表----------------------------------% 方差来源自由度平方和均方F值p值% 回归 5.0000 15.0277 3.0055 7.6122 0.0219% 残差 5.0000 1.9742 0.3948% 总计10.0000 17.0018%% 均方根误差(Root MSE) 0.6284 判定系数(R-Square) 0.8839 % 因变量均值(Dependent Mean) 4.7273 调整的判定系数(AdjR-Sq) 0.7678%% -----------------------------------参数估计-----------------------------------% 变量估计值标准误t值p值% 常数项30.9428 2011.1117 0.0154 0.9883% X1 0.7040 0.6405 1.0992 0.3218% X2 -0.8487 29.1537 -0.0291 0.9779% X1*X2 -0.0058 0.0044 -1.3132 0.2461% X1*X1 0.0003 0.0003 0.8384 0.4400% X2*X2 0.0052 0.1055 0.0492 0.9626%% Copyright 2009 - 2010 xiezhh.% $Revision: 1.0.0.0 $ $Date: 2009/12/22 21:41:00 $if nargin < 2error('至少需要两个输入参数');endp = size(X,2); % X的列数，即变量个数if nargin < 3 || isempty(model)model = 'linear'; % model参数的默认值end% 生成变量标签varnamesif nargin < 4 || isempty(varnames)varname1 = strcat({'X'},num2str([1:p]'));varnames = makevarnames(varname1,model); % 默认的变量标签elseif ischar(varnames)varname1 = cellstr(varnames);elseif iscell(varnames)varname1 = varnames(:);elseerror('varnames 必须是字符矩阵或字符串元胞数组');endif size(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,'均方根误差(Root MSE)',sqrt(ST.mse),'判定系数(R-Square)',ST.rsquare);fprintf('\n');fprintf(fmt,'因变量均值(Dependent Mean)',mean(y),'调整的判定系数(Adj R-Sq)',...ST.adjrsquare);fprintf('\n');fprintf('\n');% 显示参数估计及t检验相关结果fprintf('-------------------------------参数估计-------------------------------');fprintf('\n');fprintf('%8s%18s%15s%15s%12s','变量','估计值','标准误','t值','p值');fprintf('\n');for i = 1:size(t.beta,1)if i == 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));endendif nargout == 1stats = ST; % 返回一个包括了回归分析的所有诊断统计量的结构体变量end% -----------------子函数-----------------------function varnames = makevarnames(varname1,model)% 生成指定模型的变量标签p = size(varname1,1);varname2 = [];for i = 1:p-1varname2 = [varname2;strcat(varname1(i),'*',varname1(i+1:end))]; endvarname3 = strcat(varname1,'*',varname1);switch modelcase 'linear'varnames = varname1;case 'interaction'varnames = [varname1;varname2];case 'quadratic'varnames = [varname1;varname2;varname3];case 'purequadratic'varnames = [varname1;varname3];end调用自编的reglm函数进行二次拟合，命令如下：>> z = [9.73875 20.75 36.5987513.5725 29.6325 50.9387518.97875 36.59875 80.1387520.75125 38.22125 90.92522.055 44.58 104.7725];>> [p,t] = meshgrid([0.8 1 1.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) pp tt pp.*tt pp.^2 tt.^2]*stats.beta;>> mesh(pnew,tnew,reshape(zhat,[20,20]));>> hold on>> plot3(p(:),t(:),z(:),'k*')拟合结果：------------------------------------方差分析表------------------------------------方差来源自由度平方和均方F值p值回归 5.0000 11548.9147 2309.7829 93.4739 0.0000残差9.0000 222.3942 24.7105总计14.0000 11771.3089均方根误差(Root MSE) 4.9710 判定系数(R-Square) 0.9811因变量均值(Dependent Mean) 41.2168 调整的判定系数(AdjR-Sq) 0.9706-----------------------------------参数估计----------------------------------- 变量估计值标准误t值p值常数项242.6188 69.0439 3.5140 0.0066p -513.7781 137.3777 -3.7399 0.0046t -0.3637 0.1212 -3.0002 0.0150p*t 0.6022 0.0926 6.5010 0.0001p*p 272.2625 68.0677 3.9999 0.0031t*t -0.0003 0.0002 -1.1946 0.26282.三个自变量，一个因变量clear,clcx1=[333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 333.15 328.15 330.65 333.15 335.65 338.15 340.65 343.15 333.15 333.15 333.15 323.15 325.65 345.65 348.15]';x2=[1.19 1.206 1.228 1.23 1.252 1.27 1.277 1.31 1.35 1.39 1.43 1.23 1.23 1.23 1.23 1.23 1.2 1.2 1.2 1.2 1.2 1.26 1.26 1.26 1.26 1.26 1.23 1.23 1.23 1.23 1.23 1.23 1.23 1.15 1.47 1.51 1.23 1.23 1.23 1.23]';x3=[80 80 80 80 80 80 80 80 80 80 80 77 78 79 80 81 67 68 69 70 71 86 87 88 89 90 80 80 80 80 80 80 80 80 80 80 80 80 80 80]';y=[59.49 55.16 50.18 49.78 45.75 42.96 41.96 37.87 33.96 30.83 28.29 47.92 48.54 49.19 49.78 50.42 47.49 48.21 48.9 49.63 50.32 47.8 48.38 48.91 49.47 50.04 50.49 50.14 49.79 49.45 49.13 48.81 48.5 74.13 26.18 24.39 51.22 50.85 48.21 47.92]';X=[x1,x2,x3];ymin=min(y);y=y-ymin;fx=@(b,x1,x2,x3)(b(1)+b(2)*x1+b(3)*x2+b(4)*x3+b(5)*x1.^2+b(6)*x2.^2 +b(7)*x3.^2+b(8)*x1.*x2+b(9)*x1.*x3+b(10)*x2.*x3+b(11)*x1.^3+b(12) *x2.^3+b(13)*x3.^3)./(1+b(14)*exp(b(15)*x1+b(16)*x2+b(17)*x3+b(18) *x1.*x2+b(19)*x1.*x3+b(20)*x2.*x3));fx2=@(b,X,y)(b(1)+b(2)*X(:,1)+b(3)*X(:,2)+b(4)*X(:,3)+b(5)*X(:,1).^2+b( 6)*X(:,2).^2+b(7)*X(:,3).^2+b(8)*X(:,1).*X(:,2)+b(9)*X(:,1).*X(:,3)+b(10) *X(:,2).*X(:,3)+b(11)*X(:,1).^3+b(12)*X(:,2).^3+b(13)*X(:,3).^3)./(1+b(14)*exp(b(15)*X(:,1)+b(16)*X(:,2)+b(17)*X(:,3)+b(18)*X(:,1).*X(:,2)+b(19) *X(:,1).*X(:,3)+b(20)*X(:,2).*X(:,3)));bm=[105091.513651451,1328.10332025611,-711027.452435498,-1213.614 05762992,-1.88264106646625,934239.742471165,-25.5844409887743,-130 1.90766627356,10.5189174978167,-642.229950374061,0.002213356597694 81,-244987.606559315,0.155404373719581,9.28886223888986e-05,-0.0142 397533119651,13.4903417277274,0.0213803812532436,-0.0014125144376 6222,0.000377042917999337,-0.0845412180650883];b=bm;for l=1:10b=lsqcurvefit(fx2,b,X,y);b=nlinfit(X,y,fx2,b);endbm1=mean(x1);m2=mean(x2);m3=mean(x3);r1=range(x1); r2=range(x2);r3=range(x3);ry=range(y);x1a=min(x1);x1b=max(x1);x2a=min(x2);x2b=max(x2);x3a=min(x3);x3b=max(x3);ya=min(y);yb=max(y);n=length(y);str=num2str([1:n]');figure(1),clfplot3(x1,x2,y,'o')stem3(x1,x2,y,'filled')text(x1,x2,y+.04*ry,str,'fontsize',12)pause(.0001)hold on[x11,x22]=meshgrid(x1a:r1/75:x1b,x2a:r2/75:x2b);y1=fx(bm,x11,x22,m3);surf(x11,x22,y1)axis([x1a x1b x2a x2b ya yb])alpha(.85)shading interpaxis tightpause(1.0001)%clf% for l=1:10% plot3(x1,x2,y,'o')% stem3(x1,x2,y,'filled')% text(x1,x2,y+.04*ry,str,'fontsize',12)% pause(.0001)% hold on% m3=x3a+l*r3/10;% y1=fx(bm,x11,x22,m3);% surf(x11,x22,y1)% axis([x1a x1b x2a x2b ya yb])% alpha(.4)% shading interp% axis tight% pause(.5001)% endxlabel('X1'),ylabel('X2'),zlabel('Y')figure(2),clf[x11,x33]=meshgrid(x1a:r1/75:x1b,x3a:r3/75:x3b); plot3(x1,x3,y,'o')stem3(x1,x3,y,'filled')text(x1,x3,y+.04*ry,str,'fontsize',12)pause(.0001)hold ony2=fx(bm,x11,m2,x33);surf(x11,x33,y2)axis([x1a x1b x3a x3b ya yb])alpha(.85)shading interpaxis tightpause(5.0001)%clf% for l=1:10% plot3(x1,x3,y,'o')% stem3(x1,x3,y,'filled')% text(x1,x3,y+.04*ry,str,'fontsize',12)% pause(.0001)% hold on% m2=x2a+(l-1)*r2/10;% y2=fx(bm,x11,m2,x33);% surf(x11,x33,y2)% axis([x1a x1b x3a x3b ya yb])% alpha(.4)% shading interp% axis tight% pause(.5001)% endxlabel('X1'),ylabel('X3'),zlabel('Y')figure(3),clfplot3(x2,x3,y,'o')stem3(x2,x3,y,'filled')text(x2,x3,y+.04*ry,str,'fontsize',12)pause(.0001)hold on[x22,x33]=meshgrid(x2a:r2/75:x2b,x3a:r3/75:x3b);y3=fx(bm,m1,x22,x33);surf(x22,x33,y3)axis([x2a x2b x3a x3b ya yb])alpha(.85)shading interpaxis tightpause(5.0001)%clf% for l=1:10% plot3(x2,x3,y,'o')% stem3(x2,x3,y,'filled')% text(x2,x3,y+.04*ry,str,'fontsize',12)% pause(.0001)% hold on% m1=x1a+(l-1)*r1/10;% y3=fx(bm,m1,x22,x33);% surf(x22,x33,y3)% axis([x2a x2b x3a x3b ya yb])% alpha(.4)% shading interp% axis tight% pause(.5001)% endxlabel('X2'),ylabel('X3'),zlabel('Y')disp([' x1 x2 x3 y yhat']) yhat=fx(b,x1,x2,x3);[x1,x2,x3,y+ymin,yhat+ymin]SSy=var(y)*(n-1)RSS=(y-yhat)'*(y-yhat)Raqaure=(SSy-RSS)/SSy。