Mathematica在多元线性回归分析中的应用_黄志鹏.pdf

《Mathematica在多元线性回归分析中的应用_黄志鹏.pdf》由会员上传分享，免费在线阅读，更多相关内容在行业资料-天天文库。

1、第38卷第4期曲阜师范大学学报Vol．38No．42012年10月JournalofQufuNormalUniversityOct．2012*Mathematica在多元线性回归分析中的应用①②黄志鹏，李思泽(①北京交通大学理学院数学系;②轨道交通控制与安全国家重点实验室北京交通大学，100044，北京市)摘要:从多元线性回归分析的原理出发，利用最小二乘法准则，将回归分析归结为用Seidel迭代法求矩阵方程组解，确定待定系数的过程．利用mathematica实现精度可调节式的多元线性回归．结合实例，通过方差分

2、析表，进行F检验，计算P值，决定系数，将此方法的回归分析结果与调用统计软件包中的regress函数的结果进行比较，说明此方法的有效性．关键词:多元线性回归;最小二乘法;mathematica;Seidel迭代法中图分类号:O212文献标识码:A文章编号:1001-5337(2012)04-0028-040引言但实际背景下，可能无法找到确定的系数βi，使得上式等号成立，需要在一定准则下找到最优的待线性回归是非线性回归的基础，可以结合如神定系数，亦即找到数值解说明因变量与自变量之间［1］的线性关系．经网络的技术来

3、描述非线性系统，也可以应用将T非线性模型线性化的方法实现对非线性系统的回归假设有n组观测数据，记fi(x)=pix－bi，i=1，［2］分析，所以研究线性回归仍然是基本的和重要2，…，n，，x为待定系数βi组成的m维列向量，m也的．在线性回归方面，有学者利用mathematica研究为待定系数的个数，bi为与每组观测数据的自变量［3］了一元线性回归的mathematica实现方法，但直接组成的向量pi对应的因变量的值．最小二乘法就是使用mathematica的内置函数，不能实现编程的可调解这样的问题:节控制．

4、而且至今在多元线性回归领域，少有利用n2n+1min∑fi(x)，x∈瓗．mathematica实现可调节式编程的报道．本文旨在从i=1多元线性回归的原理出发，利用mathematica的编程用矩阵的形式表达上式，记［4］T方法，研究mathematica环境下的实现方法．ép1ùéb1ùA=êú，b=êúêúêú1回归模型系数的确定ëpTûëbûnn其中A为n×m矩阵，b为n维列向量．文献［6］给在回归分析中，如果有两个或两个以上的自变T出了当AA非奇异时，上述最小二乘问题的全局极量则称为多元回归．实际

5、中，一个因素也经常是由多_T－1T小点x=(AA)Ab即为回归模型的待定系数．个因素共同作用的结果，故多元线性回归也是一元线性回归的一般化．多元线性回归模型，为2Mathematica具体实现n+1Y=∑βiXi，出于与文献［5］中用regress函数进行对比的需i=0其中βi为待定的系数，Xi为自变量，Y为因变量，i=要，这里选取文献［5］中实验11－7的例子，进行实0，…，n+1例分析．*收稿日期:2012-03-30基金项目:轨道交通控制与安全国家重点实验室(北京交通大学)自主研究课题资助(RCS201

6、0ZT012)．作者简介:李思泽(通讯作者)，男，1957-，教授;研究方向:代数，矩阵理论，马尔科夫理论，图论在无线通信中的应用等;E-mail:lsz00@163．com．第4期黄志鹏，等:Mathematica在多元线性回归分析中的应用292．1回归模型的建立Module［2．1.1观测数据的输入{x0，eps1，eps，i，t1，j，u1，k，err，u2，nmax，x}，通过键盘输入数据组数、数据阵由m－1列的nmax=500;自变量矩阵和因变量构成的向量，得到x0=Input［"输入迭代初值向量x

7、0"］;é1x11x1(m－1)ùeps1=0.01;êúéb1ùeps=Input［"输入精度控制eps="］;A=êú，b=êúêúêúx=x0;êúëbûmë1x…xûDo［If［Abs［a［［i，i］］］＜eps1，t1=1;Return［］，t1=n1n(m－1)A中的第i行为第i次观测得到的m－1个自变量0］，{i，1，n}］;If［t1==1，Print［"Seidel迭代法失效"］，Do［Do［u1构成的行，由于要与β0这个常数项的系数配对，故=Sum［a［［i，j］］*x［［

8、j］］，{j，1，i－1}］;令A的第一列全为1(这通过对一个n维单位阵取u2=Sum［a［［i，j］］*x0［［j］］，{j，i+1，n}］;迹向量再用Join实现)．再用Transpose实现转置，最终得到ATA和ATb，即解线性方程组:ATAX=ATb．x［［i］］=(b［［i］］－u1－u2)/a［［i，i］］，{i，1，n}］;err=Max［Abs［x－x0］］;此阶段代码如下:Clea