数据统计建模方法-课件

数学建模中的统计方法 12020/12/2中国大学生数学建模竞赛（CUMCM）1992年中国工业与应用数学学会(CSIAM)开始组织 1994年起教育部高教司和CSIAM共同举办(每年9月)2010 年，来自全国33个省/市/自治区(包括香港和澳门特区)及新加坡和澳大利亚的1197所院校、17317个队（其中本科组14108队、专科组3209队）、5万多名大学生参加了本项竞赛（2001年起刊登于当年“工程数学学报”）网址：奖励：证书 （“一次参赛，终身受益”）等级：全国一等2%、二等 6%；赛区奖1/3 22020/12/2精品资料32020/12/2你怎么称呼老师？如果老师最后没有总结一节课的重点的难点，你是否会认为老师的教学方法需要改进？你所经历的课堂，是讲座式还是讨论式？教师的教鞭“不怕太阳晒，也不怕那风雨狂，只怕先生骂我笨，没有学问无颜见爹娘”“太阳当空照，花儿对我笑，小鸟说早早早”42020/12/2我国CUMCM竞赛规模52020/12/2学生欢迎：“一次参赛，终身受益”研究生导师们的认同企业界的认同赞助教育改革同行的认同：“成功范例”国际同行的认同竞赛的反响62020/12/2IBM 中国研究中心-招聘条件Position title:Business Optimization(BJ)1Background in industrial engineering,operations research,mathematics,Artificial Intelligence,management science etc.2.Knowledge in network design,job scheduling,data analysis,simulation and optimization 3.Award in mathematical contest in modeling is a plus 4.Experience in industry is a plus 5.Experience in eclipse or programming model/architecture design is a plus-Feb.18,2006,竞赛的反响（一例）72020/12/2IBM 中国研究中心:Business Analysis Optimization Job Requirements:1、PhD M.S.in mathematics,statistics,computer science,industrial engineering management science etc.2、Self-motivated,responsible,able to wk independently under tight deadline willing to wk under pressure.3、Skill in applied mathematics,including mathematical programming,statistics,data mining,simulation etc.4、Knowledge in supply chain logistics strategy modeling,simulation,planning optimization.5、Strong interest basic knowledge about industry trends,technologies,solutions in analytics optimization.6、Experience in ERP/SCM/CRM system SCM consulting practice is a plus.7、Award in highly regarded mathematical modeling contest is a plus.8、Experience in eclipse,Java,architecture design is a plus.-March 26,2009,竞赛的反响（一例）82020/12/2CUMCM评阅标准清晰性：摘要应理解为详细摘要，提纲挈领 表达严谨、简捷，思路清新 格式符合规范，严禁暴露身份创造性：特别欣赏独树一帜、标新立异，但要合理假设的合理性，建模的创造性，结果的正确性，表述的清晰性。正确性：不强调与“参考答案”的一致性和结果的精度；好方法的结果一般比较好；但不一定是最好的合理性：关键假设(不欣赏罗列大量无关紧要的假设);要对假设的合理性进行解释，正文中引用92020/12/2CUMCM评阅标准:一些常见问题有的论文过于简单，该交代的内容省略了，难以看懂有的队罗列一系列假设或模型，又不作比较、评价，希望碰上“参考答案”或“评阅思路”，弄巧成拙数学模型最好明确、合理、简洁：有些论文不给出明确的模型，只是根据赛题的情况，实际上是用“凑”的方法给出结果，虽然结果大致是对的，没有一般性，不是数学建模的正确思路。有的论文参考文献不全，或引用他人结果不作交代;参考文献应在正文中引用102020/12/2从论文评阅看学生参加竞赛中的问题 吃透题意方面不足，没有抓住和解决主要问题；就事论事，形成数学模型的意识和能力欠缺；对所用方法一知半解，不管具体条件，套用现成的方法，导致错误；对结果的分析不够，怎样符合实际考虑不周；写作方面的问题(摘要、简明、优缺点、参考文献);队员之间合作精神差，孤军奋战；依赖心理重，甚至违纪（指导教师、网络）。112020/12/2竞赛内容与形式内容 赛题：工程、管理中经过简化的实际问题 答卷：一篇包含问题分析、模型假设、建立、求解(通常用计算机)、结果分析和检验等的论文形式 3名大学生组队，在3天内完成的通讯比赛 可使用任何“死”材料(图书/互联网/软件等),但不得与队外任何人讨论（包括上网讨论）宗旨创新意识 团队精神 重在参与 公平竞争标准假设的合理性，建模的创造性，结果的正确性，表述的清晰性。122020/12/22020/12/213一、一、CUMCM历年赛题的分析历年赛题的分析3 3、从问题的解决方法上分析、从问题的解决方法上分析 从问题的解决方法上分析，涉及到的数学建从问题的解决方法上分析，涉及到的数学建模方法：模方法：几何理论、组合概率、几何理论、组合概率、统计统计(回归回归)分分析析、优化方法（规划）、优化方法（规划）、图论与网络优化图论与网络优化、层、层次分析、插值与拟合、差分方法、微分方程、次分析、插值与拟合、差分方法、微分方程、排队论、模糊数学、随机决策、多目标决策、排队论、模糊数学、随机决策、多目标决策、随机模拟、灰色系统理论、神经网络、时间序随机模拟、灰色系统理论、神经网络、时间序列、综合评价、机理分析等方法。列、综合评价、机理分析等方法。2020/12/214 用的最多的方法是优化方法和概率统计的方法用的最多的方法是优化方法和概率统计的方法.用到优化方法的共有用到优化方法的共有2121个题，占总数的个题，占总数的70%70%，其，其中整数规划中整数规划4 4个，线性规划个，线性规划6 6个，非线性规划个，非线性规划1414个个,多目标规划多目标规划5 5个。个。用到概率统计方法的有用到概率统计方法的有1616个题，占个题，占53.3%53.3%，几乎，几乎是每年至少有一个题目用到概率统计的方法。是每年至少有一个题目用到概率统计的方法。用到图论与网络优化方法的问题有用到图论与网络优化方法的问题有5 5个；个；用到层次分析方法的问题有个；用到层次分析方法的问题有个；一、一、CUMCM历年赛题的分析历年赛题的分析3 3、从问题的解决方法上分析、从问题的解决方法上分析2020/12/215 用到插值拟合的问题有用到插值拟合的问题有6 6个；个；用到神经网络的用到神经网络的4 4个；个；用灰色系统理论的用灰色系统理论的4 4个个;用到时间序列分析的至少用到时间序列分析的至少2 2个个;用到综合评价方法的至少用到综合评价方法的至少3 3个；个；机理分析方法和随机模拟都多次用到机理分析方法和随机模拟都多次用到;其他的方法都至少用到一次。其他的方法都至少用到一次。大部分题目都可以用两种以上的方法来解决大部分题目都可以用两种以上的方法来解决,即综合性较强的题目有即综合性较强的题目有2525个，占个，占83.3%83.3%。一、一、CUMCM历年赛题的分析历年赛题的分析3 3、从问题的解决方法上分析、从问题的解决方法上分析统计的基本概念统计的基本概念参数估计参数估计假设检验假设检验数数据据的的统统计计描描述述和和分分析析162020/12/2一、统计量一、统计量172020/12/2182020/12/2二、分布函数的近似求法二、分布函数的近似求法192020/12/2三、几个在统计中常用的概率分布三、几个在统计中常用的概率分布-4-2024600.050.10.150.20.250.30.350.41正态分布正态分布),(2smN密度函数：222)(21)(smsp-=xexp分布函数：dyexFyx222)(21)(smsp-=其中m为均值，2s为方差，+-x.标准正态分布：N（0，1）密度函数2221)(xex-=pjdyexyx2221)(-=Fp 分布函数202020/12/2212020/12/2222020/12/2232020/12/2242020/12/2252020/12/2F（10，50）分布的密度函数曲线262020/12/2统计工具箱中的基本统计命令统计工具箱中的基本统计命令1.数据的录入、保存和调用数据的录入、保存和调用2.基本统计量基本统计量3.常见的概率分布函数常见的概率分布函数4.频频 数数 直直 方方 图图 的的 描描 绘绘5.参数估计参数估计6.假设检验假设检验7.综合实例综合实例272020/12/2一、数据的录入、保存和调用一、数据的录入、保存和调用 例例1 上海市区社会商品零售总额和全民所有制职工工资总额的数据如下：统计工具箱中的基本统计命令统计工具箱中的基本统计命令282020/12/21年份数据以1为增量，用产生向量的方法输入.命令格式：x=a:h:b t=78:872分别以x和y代表变量职工工资总额和商品零售总额.x=23.8,27.6,31.6,32.4,33.7,34.9,43.2,52.8,63.8,73.4 y=41.4,51.8,61.7,67.9,68.7,77.5,95.9,137.4,155.0,175.03将变量t、x、y的数据保存在文件data中.save data t x y 4进行统计分析时，调用数据文件data中的数据.load data292020/12/21输入矩阵：data=78,79,80,81,82,83,84,85,86,87,88;23.8,27.6,31.6,32.4,33.7,34.9,43.2,52.8,63.8,73.4;41.4,51.8,61.7,67.9,68.7,77.5,95.9,137.4,155.0,175.02将矩阵data的数据保存在文件data1中：save data1 data3 3进行统计分析时，先用命令：load data1 调用数据文件data1中的数据，再用以下命令分别将矩阵data的第一、二、三行的数据赋给变量t、x、y：t=data(1,:)x=data(2,:)y=data(3,:)若要调用矩阵data的第j列的数据，可用命令：data(:,j)302020/12/2二、基本统计量二、基本统计量对随机变量x，计算其基本统计量的命令如下：均值：mean(x)中位数：median(x)标准差：std(x)方差：var(x)偏度：skewness(x)峰度：kurtosis(x)例例 对例1中的职工工资总额x，可计算上述基本统计量.312020/12/2 例例1 上海市区社会商品零售总额和全民所有制职工工资总额的数据如下：x=23.8 27.6 31.6 32.4 33.7 34.9 43.2 52.8 63.8 73.4;y=41.4 51.8 61.7 67.9 68.7 77.5 95.9 137.4 155.0 175.0;X=mean(x),median(x),std(x),var(x),skewness(x),kurtosis(x)X=41.7200 34.3000 16.4544 270.7462 0.8590 2.4037Y=mean(y),median(y),std(y),var(y),skewness(y),kurtosis(y)Y=1.0e+003*0.0932 0.0731 0.0464 2.1497 0.0007 0.0020322020/12/2三三、常见概率分布的函数常见概率分布的函数MATLAB工具箱对每一种分布都提供5类函数，其命令字符为：概率密度：pdf 概率分布：cdf逆概率分布：inv 均值与方差：stat随机数生成：rnd （当需要一种分布的某一类函数时，将以上所列的分布命令字符与函数命令字符接起来，并输入自变量（可以是标量、数组或矩阵）和参数即可.）332020/12/2在MATLAB中输入以下命令：x=-6:0.01:6;y=normpdf(x);z=normpdf(x,0,2);plot(x,y,x,z)1密度函数密度函数：p=normpdf(x,mu,sigma)(当mu=0,sigma=1时可缺省)如对均值为mu、标准差为sigma的正态分布，举例如下：342020/12/23逆概率分布逆概率分布：x=norminv(P,mu,sigma).即求出x，使得PX50），按中心极限定理，它近似地 服从正态分布；二、使用MATLAB工具箱中具有特定分布总体的估计命令.（1）muhat,muci=expfit(X,alpha)在显著性水平alpha下，求指数分布的数据X的均值的点估计及其区间估计.（2）lambdahat,lambdaci=poissfit(X,alpha)在显著性水平alpha下，求泊松分布的数据X的参数的点估计及其区间估计.（3）phat,pci=wblfit(X,alpha)在显著性水平alpha下，求Weibull分布的数据X的参数的点估计及其区间估计.402020/12/2六、假设检验六、假设检验 在总体服从正态分布的情况下，可用以下命令进行假设检验.1总体方差总体方差 已知时，总体均值的检验使用已知时，总体均值的检验使用 z检验检验 h,sig,ci=ztest(x,m,sigma,alpha,tail)检验数据 x 的关于均值的某一假设是否成立，其中sigma 为已知方差，alpha 为显著性水平，究竟检验什么假设取决于 tail 的取值：tail=0，检验假设“x 的均值等于等于 m”tail=1，检验假设“x 的均值大于大于 m”tail=-1，检验假设“x 的均值小于小于 m”tail的缺省值为 0，alpha的缺省值为 0.05.返回值 h 为一个布尔值，h=1 表示可以拒绝假设，h=0 表示不可以拒绝假设，sig 为假设成立的概率，ci 为均值的 1-alpha 置信区间.412020/12/2 例例 MATLAB统计工具箱中的数据文件gas.mat.中提供了美国1993年1月份和2月份的汽油平均价格（price1,price2分别是1、2月份的油价，单位为美分），它是容量为20的双样本.假设1月份油价的标准偏差是每加仑4分币（=4），试检验1月份油价的均值是否等于115.解解 作假设：m=115.首先取出数据，用以下命令：load gas然后用以下命令检验 h,sig,ci=ztest(price1,115,4)返回：h=0，sig=0.8668，ci=113.3970 116.9030.检验结果:1.布尔变量h=0,表示不拒绝零假设.说明提出的假设均值115 是合理的.2.sig值为0.8668,远超过0.5,不能拒绝零假设 3.95%的置信区间为113.4,116.9,它完全包括115,且精度很 高.422020/12/22总体方差总体方差 未知时，总体均值的检验使用未知时，总体均值的检验使用t 检验检验 h,sig,ci=ttest(x,m,alpha,tail)检验数据 x 的关于均值的某一假设是否成立，其中alpha 为显著性水平，究竟检验什么假设取决于 tail 的取值：tail=0，检验假设“x 的均值等于等于 m”tail=1，检验假设“x 的均值大于大于 m”tail=-1，检验假设“x 的均值小于小于 m”tail的缺省值为 0，alpha的缺省值为 0.05.返回值 h 为一个布尔值，h=1 表示可以拒绝假设，h=0 表示不可以拒绝假设，sig 为假设成立的概率，ci 为均值的 1-alpha 置信区间.432020/12/2返回：h=1，sig=4.9517e-004，ci=116.8 120.2.检验结果:1.布尔变量h=1,表示拒绝零假设.说明提出的假 设油价均值115是不合理的.2.95%的置信区间为116.8 120.2,它不包括 115,故不能接受假设.3.sig值为4.9517e-004,远小于0.5,不能接受零 假设.例例8 试检验例8中2月份油价price2的均值是否等于115.解解 作假设：m=115，price2为2月份的油价，不知其方差，故用以下命令检验h,sig,ci=ttest(price2,115)442020/12/23两总体均值的假设检验两总体均值的假设检验使用使用 t 检验检验 h,sig,ci=ttest2(x,y,alpha,tail)检验数据 x，y 的关于均值的某一假设是否成立，其中alpha 为显著性水平，究竟检验什么假设取决于 tail 的取值：tail=0，检验假设“x 的均值等于等于 m”tail=1，检验假设“x 的均值大于大于 m”tail=-1，检验假设“x 的均值小于小于 m”tail的缺省值为 0，alpha的缺省值为 0.05.返回值 h 为一个布尔值，h=1 表示可以拒绝假设，h=0 表示不可以拒绝假设，sig 为假设成立的概率，ci 为与x与y均值差的的 1-alpha 置信区间.452020/12/2返回：h=1，sig=0.0083，ci=-5.8,-0.9.检验结果:1.布尔变量h=1,表示拒绝零假设.说明提出的 假设“油价均值相同”是不合理的.2.95%的置信区间为-5.8,-0.9,说明一月份油 价比二月份油价约低1至6分.3.sig-值为0.0083,远小于0.5,不能接受“油价均 相同”假设.例例9 试检验例8中1月份油价price1与2月份的油价price2均值是否相同.解解 用以下命令检验h,sig,ci=ttest2(price1,price2)462020/12/24非参数检验：总体分布的检验非参数检验：总体分布的检验MATLAB工具箱提供了两个对总体分布进行检验的命令:（1）h=normplot(x)（2）h=wblplot(x)此命令显示数据矩阵x的正态概率图.如果数据来自于正态分布，则图形显示出直线性形态.而其它概率分布函数显示出曲线形态.此命令显示数据矩阵x的Weibull概率图.如果数据来自于Weibull分布，则图形将显示出直线性形态.而其它概率分布函数将显示出曲线形态.472020/12/2例例10 一道工序用自动化车床连续加工某种零件，由于刀具损坏等会出现故障.故障是完全随机的，并假定生产任一零件时出现故障机会均相同.工作人员是通过检查零件来确定工序是否出现故障的.现积累有100次故障纪录，故障出现时该刀具完成的零件数如下：x=459 362 624 542 509 584 433 748 815 505.612 452 434 982 640 742 565 706 593 680.926 653 164 487 734 608 428 1153 593 844.527 552 513 781 474 388 824 538 862 659.402 960 885 610 292 837 473 677 358 638.699 634 555 570 84 416 606 1062 484 120.447 654 564 339 280 246 687 539 790 581.621 724 531 512 577 496 468 499 544 645.764 558 378 765 666 763 217 715 310 851;试观察该刀具出现故障时完成的零件数属于哪种分布.482020/12/2解解 1数据输入2作频数直方图 hist(x,10)3分布的正态性检验 normplot(x)4参数估计：muhat,sigmahat,muci,sigmaci=normfit(x)（看起来刀具寿命服从正态分布）估计出该刀具的均值为587，方差199，均值的0.95置信区间为544.9090，628.1577，方差的0.95置信区间为173.3397，232.9193.492020/12/25假设检验 已知刀具的寿命服从正态分布，现在方差未知的情况下，检验其均值 m 是否等于587.结果：h=0，sig=0.9823，ci=544.9090，628.1577.检验结果:1.布尔变量h=0,表示不拒绝零假设.说 明提出的假设寿命均值587是合理的.2.95%的置信区间为544.9，628.16,它完全包括587,且精度很高.3.sig值为0.9823,远超过0.5,不能拒绝 零假 设.502020/12/2一元线性回归一元线性回归多元线性回归多元线性回归二、回归分析二、回归分析方法方法数数学学模模型型及及定定义义模模型型参参数数估估计计检检验验、预预测测与与控控制制可可线线性性化化的的一一元元非非线线性性回回归归（曲曲线线回回归归）数数学学模模型型及及定定义义模模型型参参数数估估计计多多元元线线性性回回归归中中的的检检验验与与预预测测512020/12/2522020/12/22.12.1一一元线性回归分析元线性回归分析例例1 测16名成年女子的身高与腿长所得数据如下：以身高x为横坐标，以腿长y为纵坐标将这些数据点（xi，yi）在平面直角坐标系上标出.散点图身高身高143145146147149150153154155156157158159160162164腿长腿长88858891929393959698979698991001022.1.1 数学模型数学模型532020/12/2一元线性回归分析的主要任务是：一元线性回归分析的主要任务是：、1、用试验值、用试验值（样本值）对（样本值）对0b b1b b和和s s作点估计；作点估计；、2、对回归系数、对回归系数0b b1b b作假设检验；作假设检验；542020/12/22.2.2 模型参数估计模型参数估计（1）回归系数的最小二乘估计）回归系数的最小二乘估计552020/12/2其中其中 562020/12/2572020/12/2582020/12/2称Qe为残差平方和残差平方和或剩余平方和剩余平方和.可以证明：592020/12/22.2.3 2.2.3 回归方程的显著性检验回归方程的显著性检验602020/12/2612020/12/2且 与 独立1t检验法检验法若H0成立，即b=0，可以用三种不同方法进行检验，它们的本质可以用三种不同方法进行检验，它们的本质是相同的这里介绍是相同的这里介绍622020/12/2因而故为显著性水平即得H0的拒绝域为632020/12/2642020/12/2 2.2.F F检验法检验法：当当H H0 0为真时，为真时，且且SSR与与SSE相互独立；因此，当相互独立；因此，当H0为真时，为真时，当当FFFF1-1-(1,n-2)(1,n-2)时应该放弃原假设时应该放弃原假设H H0 0。652020/12/23.r-检验法662020/12/2回归系数的置信区间回归系数的置信区间672020/12/2线性回归线性回归matlabmatlab b=regress(Y,X)1确定回归系数的点估计值：确定回归系数的点估计值：b为一元函数为一元函数y系数的估计值系数的估计值Y.X都是列向量682020/12/22求回归系数的点估计和区间估计、并检验回归模型：求回归系数的点估计和区间估计、并检验回归模型：b,bint,r,rint,stats=regress(Y,X,alpha)回归系数的区间估计残差用于检验回归模型的统计量，有三个数值：相关系数R 2、F值、与F 对应的概率p置信区间 显著性水平（缺省时为0.05）692020/12/23画出残差及其置信区间：画出残差及其置信区间：rcoplot（r，rint）702020/12/2例例1 测16名成年女子的身高与腿长所得数据如下：身高身高143145146147149150153154155156157158159160162164腿长腿长8885889192939395969897969899100102clearx=143 145 146 147 149 150 153 154 155 156 157 158 159 160 162 164;X=ones(16,1)x;Y=88 85 88 91 92 93 93 95 96 98 97 96 98 99 100 102;b,bint,r,rint,stats=regress(Y,X);b,bint,statsrcoplot(r,rint)z=b(1)+b(2)*x;figure,plot(x,Y,k+,x,z,r)F=finv(0.95,1,14)%及求F0.05(1,16-2)712020/12/2b=-16.0730 0.7194bint=-33.7071 1.5612 0.6047 0.8340stats=0.9282 180.9531 0.0000F0.05(1,16-2)=4.6722020/12/2 从残差图可以看出，除第二个数据外，其余数据的残差离零点均较近，且残差的置信区间均包含零点，这说明回归模型 y=-16.073+0.7194x能较好的符合原始数据，而第二个数据可视为异常点.732020/12/22.2.4 预测与控制预测与控制（1）预测）预测:对固定的对固定的x值预测相应的值预测相应的y值值742020/12/2例例1 测16名成年女子的身高与腿长所得数据如下：身高身高143145146147149150153154155156157158159160162164腿长腿长8885889192939395969897969899100102请你预测一下身高为请你预测一下身高为170cm时该女子的腿长应该是在什么范围时该女子的腿长应该是在什么范围(置信度为置信度为0.05)clearx0=170 x=143 145 146 147 149 150 153 154 155 156 157 158 159 160 162 164;y0=-16.073+0.7194*x0;t=tinv(0.975,16-2);Lxx=sum(x-mean(x).2)d=sum(r.2)/14*t*sqrt(1+1/16+(x0-mean(x)/Lxx)yy=y0-d,y0+d752020/12/2clearx=143 145 146 147 149 150 153 154 155 156 157 158 159 160 162 164;X=ones(16,1)x;Y=88 85 88 91 92 93 93 95 96 98 97 96 98 99 100 102;b,bint,r,rint,stats=regress(Y,X);x0=170 x=143 145 146 147 149 150 153 154 155 156 157 158 159 160 162 164;y0=-16.073+0.7194*x0;t=tinv(0.975,16-2);Lxx=sum(x-mean(x).2)d=sum(r.2)/14*t*sqrt(1+1/16+(x0-mean(x)/Lxx)yy=y0-d,y0+dyy=102.3216 110.1284762020/12/2（2）控制）控制772020/12/2782020/12/2用一条曲线来拟合因变量对于自变量的依赖关系用一条曲线来拟合因变量对于自变量的依赖关系通过变量替换把问题转化为一元或多元线性回归问通过变量替换把问题转化为一元或多元线性回归问题后，用线性回归分析的方法建立回归模型，并进题后，用线性回归分析的方法建立回归模型，并进行预测行预测(即化非线性回归为线性回归即化非线性回归为线性回归)2.2 可线性化的一元非线性回归曲线回归可线性化的一元非线性回归曲线回归792020/12/21.用用幂函数曲线拟合幂函数曲线拟合 设则及化非线性回归为线性回归的常用方法802020/12/22.指数函数曲线拟合指数函数曲线拟合设则812020/12/23.3.对数函数曲线拟合对数函数曲线拟合 设则822020/12/24.4.双曲线函数拟合双曲线函数拟合 设 则832020/12/25.5.二次多项式及三次多项式二次多项式及三次多项式 842020/12/26、S型（Logistic）曲线拟合 变形 设852020/12/2例例2 出钢时所用的盛钢水的钢包，由于钢水对耐火材料的侵蚀，容积不断增大.我们希望知道使用次数与增大的容积之间的关 系.对一钢包作试验，测得的数据列于下表：862020/12/2散点图此即非线性回归非线性回归或曲线回归曲线回归问题（需要配曲线）配曲线的一般方法是：配曲线的一般方法是：872020/12/2通常选择的六类曲线如下：882020/12/2clear all,clc;x=2:1:16;x1=1./x;x1=x1;X=ones(15,1),x1;y=6.42 8.20 9.58 9.50 9.70 10.00 9.93 9.99 10.49 10.59 10.60 10.80 10.60 10.90 10.76;y1=log(y);b,bint,r,rin,stats=regress(y1,X);b,bint,statsb=2.4578 -1.1107bint=2.4306 2.4850 -1.2485 -0.9729stats=0.9589 303.1896 0.0000y=a 892020/12/22.3 多元线性回归多元线性回归多元线性回归在工程上更为有用。2.3.12.3.1数学模型及定义数学模型及定义902020/12/2912020/12/2922020/12/2932020/12/22.3.22.3.2模型参数估计模型参数估计 942020/12/2解得估计值952020/12/2962020/12/22.3.4多元线性回归中的检验与预测多元线性回归中的检验与预测972020/12/2982020/12/2992020/12/21002020/12/22.4.1多元线性回归多元线性回归 b=regress(Y,X)1)确定回归系数的点估计值：确定回归系数的点估计值：2.4 MATLAB统计工具箱中的回归分析命令对一元线性回归，取p=1即可.1012020/12/23、画出残差及其置信区间：画出残差及其置信区间：rcoplot（r，rint）2)求回归系数的点估计和区间估计、并检验回归模型：求回归系数的点估计和区间估计、并检验回归模型：b,bint,r,rint,stats=regress(Y,X,alpha)回归系数的区间估计残差用于检验回归模型的统计量，有三个数值：相关系数r2、F值、与F对应的概率p置信区间 显著性水平（缺省时为0.05）1022020/12/22.4.2多多 项项 式式 回回 归归（1）一元多项式回归）一元多项式回归 1）确定多项式系数的命令：p,S=polyfit（x，y，m）A、回归：、回归：y=a1xm+a2xm-1+amx+am+11032020/12/2B、预测和预测误差估计：、预测和预测误差估计：（1）Y=polyval（p，x）求polyfit所得的回归多项式在x处的预 测值Y；（2）Y，DELTA=polyconf（p，x，S，alpha）求polyfit所得 的回归多项式在x处的预测值Y及预测值的显著性为1-alpha的置信区间Y DELTA；alpha缺省时为0.05.一元多项式回归也可以化为多元线性回归来解。（2）一元多项式回归命令：polytool（x，y，m）此命令产生一个交互式的画面，画面中有拟合曲线和此命令产生一个交互式的画面，画面中有拟合曲线和y y的置信区间。通过左下方的的置信区间。通过左下方的ExportExport菜单，可以输出回归系数菜单，可以输出回归系数等。等。1042020/12/2法一法一 直接作二次多项式回归：直接作二次多项式回归：t=1/30:1/30:14/30;s=11.86 15.67 20.60 26.69 33.71 41.93 51.13 61.49 72.90 85.44 99.08 113.77 129.54 146.48;p,S=polyfit(t,s,2)得回归模型为：Y,DELTA=polyconf(p,t,S);plot(t,s,k+,t,Y,r)预测及作图预测及作图1052020/12/2法二法二化为多元线性回归：化为多元线性回归：t=1/30:1/30:14/30;s=11.86 15.67 20.60 26.69 33.71 41.93 51.13 61.49 72.90 85.44 99.08 113.77 129.54 146.48;T=ones(14,1)t(t.2);b,bint,r,rint,stats=regress(s,T);b,stats得回归模型为：1062020/12/2polytool(x,s,2)法三法三1072020/12/22.4.3多元二项式回归多元二项式回归命令：rstool（x，y，model,alpha）nm矩阵显著性水平（缺省时为0.05）n维列向量1082020/12/2 命令rstool产生一个交互式画面，画面中有m个图形，这m个图形分别给出了一个独立变量xi（另m-1个变量取固定值）与y的拟合曲线，以及y的置信区间。可以通过键入不同的xi值来获得相应的y值。1092020/12/2 例例3 设某商品的需求量与消费者的平均收入、商品价格的统计数 据如下，建立回归模型，预测平均收入为800、价格为6时 的商品需求量.解解 直接用多元二项式回归：x1=1000 600 1200 500 300 400 1300 1100 1300 300;x2=5 7 6 6 8 7 5 4 3 9;y=100 75 80 70 50 65 90 100 110 60;x=x1 x2;rstool(x,y,purequadratic)1102020/12/2 在画面左下方的下拉式菜单中选”all”,则beta（回归系数）、rmse（剩余标准差）和residuals（残差）都传送到Matlab工作区中.在左边图形下方的方框中输入800，右边图形下方的方框中输入6。则画面左边的“Predicted Y”下方的数据变为86.3971，即预测出平均收入为800、价格为6时的商品需求量为86.3971.1112020/12/2在Matlab工作区中输入命令：beta,rmse1122020/12/22.5非线性非线性回回 归归（1）确定回归系数的命令：beta，r，J=nlinfit（x，y，model,beta0）（2）非线性回归命令：nlintool（x，y，model,beta0，alpha）2.5.1回归：回归：残差Jacobian矩阵，用于估计预测误差需要的数据。回归系数的初值是事先用m-文件定义的非线性函数估计出的回归系数输入数据x、y分别为 矩阵和n维列向量，对一元非线性回归，x为n维列向量。其中个参数含义同前，alpha为显著性水平，缺省时为0.05。该命令产生一个交互式的画面，画面中有拟合曲线和y的置信区间。通过左下方的Export菜单，可以输出回归系数等。1132020/12/22.5.2、预测和预测误差估计：、预测和预测误差估计：该命令用于求nlinfit 或nlintool所得的回归函数在x处的预测值Y及预测值的显著性为1-alpha的置信区间Y DELTA.Y，DELTA=nlpredci（model,x，beta，r，J）1142020/12/2例例 4 对第一节例2，求解如下：2、输入数据：x=2:16;y=6.42 8.20 9.58 9.5 9.7 10 9.93 9.99 10.49 10.59 10.60 10.80 10.60.10.90 10.76;beta0=8 2;3、求回归系数：beta,r,J=nlinfit(x,y,volum,beta0);beta得结果：beta=11.6036 -1.0641即得回归模型为：题目1152020/12/24、预测及作图：YY,delta=nlpredci(volum,x,beta,r,J)；plot(x,y,k+,x,YY,r)1162020/12/2逐步回归分析逐步回归分析（Stepwise regresssion analysis）1172020/12/21182020/12/21192020/12/21202020/12/21212020/12/21222020/12/21232020/12/21242020/12/21252020/12/21262020/12/21272020/12/21282020/12/21292020/12/21302020/12/21312020/12/21322020/12/21332020/12/21342020/12/21352020/12/21362020/12/21372020/12/2逐步回归逐步回归逐步回归的命令是：stepwise(x,y,inmodel，penter,premove）运行stepwise命令时产生三个图形窗口：Stepwise Plot，Stepwise Table，Stepwise History.在Stepwise Plot窗口，显示出各项的回归系数.Stepwise Table 窗口中列出了一个统计表，包括回归系数，以及模型的统计量剩余标准差（RMSE）、相关系数（R-square）、F值、与F对应的概率P.矩阵的列数的指标，给出初始模型中包括的子集（缺省时设定为全部自变量）进入容忍值自变量数据,阶矩阵因变量数据，阶矩阵删除容忍值1382020/12/2例例6 6 水泥凝固时放出的热量y与水泥中4种化学成分x1、x2、x3、x4有关，今测得一组数据如下，试用逐步回归法确定一个线性模型.1 1、数据输入：、数据输入：x1=7 1 11 11 7 11 3 1 2 21 1 11 10;x2=26 29 56 31 52 55 71 31 54 47 40 66 68;x3=6 15 8 8 6 9 17 22 18 4 23 9 8;x4=60 52 20 47 33 22 6 44 22 26 34 12 12;y=78.5 74.3 104.3 87.6 95.9 109.2 102.7 72.5 93.1 115.9.83.8 113.3 109.4;x=x1 x2 x3 x4;stepwise(x,y,4,0.06,0.1)1392020/12/2注意注意1402020/12/21412020/12/21422020/12/21432020/12/21442020/12/2四四 软件开发人员的薪金软件开发人员的薪金问题：问题：一家高技术公司人事部门为研究软件开发人员的薪金与他们的资历、管理责任、教育程度等因素之间的关系，要建立一个数学模型，以便分析公司人士策略的合理性，并作为新聘用人员工资的参考。他们认为目前公司人员的薪金总体上是合理的，可以作为建模的依据，于是调查了46名开发人员的档案资料，如表。其中资历一列指从事专业工作的年数，管理一列中1表示管理人员，0表示非管理人员，教育一列中1表示中学程度，2表示大学程度，3表示更高程度（研究生）1452020/12/2编号薪金资历管理教育编号薪金资历管理教育0113876111131980031302116081031411417401031870111315202634130411283102161323140305117671031712844402062087221218132455020711772202191367750308105352012015965511091219520321123666011012313302222135261311149753112313839602122137131224228846121462020/12/2编号薪金资历管理教育编号薪金资历管理教育251697871136168821202261480380237241701213271740481138159901301282218481339263301312291354880140179491402301446710014125685151331159421002422783716123223174101343188381602332378010124417483160134254101112451920717023514861110146193462001开发人员的薪金与他们的资历、管理责任、教育程度1472020/12/2分析与假设：分析与假设：按照常识，薪金自然按照资历（年）的增长而增加，管理人员的薪金高于非管理人员，教育程度越高薪金越高。薪金记作，资历（年）记作，为了表示是否为管理人员定义1，管理人员0，非管理人员为了表示三种教育程度，定义1，中学0，其它1，大学0，其它这样，中学用表示，大学用表示，研究生则用表示。1482020/12/2为了简单起见，我们假定资历（年）对薪金的作用是线性的，即资历每加一年，薪金的增长是常数；管理责任、教育程度、资历诸因素之间没有交互作用，建立线性回归模型。基本模型：基本模型：薪金与资历，管理责任，教育程度之间的多元线性回归模型为其中，是待估计的回归系数，是随机误差。利用MATLAB的系统工具箱可以得到回归系数及其置信区间（置信水平 ）、检验统计量的结果，见表。1492020/12/2模型求解模型求解参数参数参数估计值参数估计值置信区间置信区间a011032 10258 11807 a1546 484 608 a26883 6248 7517 a3-2994-3826 -2162 a4148-636 931 R2=0.957 F=226 p=0.000R2,F,p 模型整体上模型整体上可用可用,Finv(0.95,4,46-4-1)=2.6资历增加资历增加1年薪金增年薪金增长长546 管理人员薪金多管理人员薪金多6883 中学程度薪金比更高的少中学程度薪金比更高的少2994 大学程度薪金比更高的多大学程度薪金比更高的多148 a4置信区间包含零点，解释置信区间包含零点，解释不可靠不可靠!中学：中学：x3=1,x4=0;大学：大学：x3=0,x4=1;更高：更高：x3=0,x4=0.x2=1 管理，管理，x2=0 非非管理管理x1资历资历(年年)注意：上述解释是就平均值来说的，并且，一个因素改变引起的因变量的变化量，都是在其它因素不变的条件下才成立的。1502020/12/21512020/12/2进一步讨论：进一步讨论：的置信区间包含零点，说明上述基本模型存在缺点。为了寻找改进的方向，常用残差分析法（残差指薪金的实际值与模型估计的薪金之差，是基本模型中随机误差的估计值，这里用同一个符号）。我们将影响因素分成资历与管理教育组合两类，管理-教育组合定义如表。组合1 23456管理010101教育112233管理管理教育组合教育组合1522020/12/2为了对残差进行分析，下图给出与资历的关系，及与管理-教育组合间的关系。从左图看，残差大概分成3个水平，这是由于6种管理教育组合混在一起，在模型中未被正确反映的结果；从右图看，对于前4个管理教育组合，残差或者全为正，或者全为负，也表明管理-教育组合在模型中处理不当。在模型中，管理责任和教育程度是分别起作用的，事实上，二者可能起着交互作用，如大学程度的管理人员的薪金会比二者分别的薪金之和高一点。e 与资历与资历x1的关系的关系 e与管理与管理教育组合的关系教育组合的关系 1532020/12/2以上分析提示我们，应在基本模型中增加管理更好的模型：更好的模型：与教育的交互项，建立新的回归模型。增加与的交互项后，模型记作利用MATLAB的统计工具箱regress得到的结果如表：1542020/12/2参数参数估计值置信区间1120411044 11363497486 50870486841 7255-1727-1939 -1514-348-545 -152-3071-3372 -276918361571 21011552020/12/21562020/12/2由上表可知，这个模型的做该模型的两个残差分析图，可以看出，已经消除了不正常现象，这也说明了模型的适用性。和值都比上一个模型有所改进，并且所有回归系数的置信区间都不含零点，表明这个模型完全可用。与的关系与组合的关系1572020/12/2从上图，还可以发现一个异常点：具有10年资历、大学程度的管理人员（编号33）的实际薪金明显低于模型的估计值，也明显低于与他有类似经历的其他人的薪金。这可能是由我们未知的原因造成的。为了使个别数据不致影响整个模型，应该将这个异常数据去掉，对模型重新估计回归系数，得到的结果如表。残差分析见图。可以看到，去掉异常数据后结果又有改善。1582020/12/2参数参数估计值置信区间1120011139 11261498494 50370416962 7120-1737-1818 -1656-356-431 -281-3056-3171 -294219971894 21001592020/12/2与的关系与组合的关系e x1 e 组合组合1602020/12/2评注：评注：从建立回归模型的角度，通过这个问题的求解我们学习了：1）对于影响因变量的定性因素（管理、教育），可以引入 01变量来处理，01变量的个数比定性因素的水平少 1（如教育程度有3个水平，引入2个01变量）。2）用残差分析法可以发现模型的缺陷，引入交互作用项常 常可以得到改善。3）若发现异常值应剔除，有助于结果的合理性。思考：思考：在这里我们由简到繁，先分别引进管理和教育因素，再引入交互项。试直接对6种管理-教育组合引入5个01变量，建立模型，看结果如何。1612020/12/2回归模型方法：大样本的内部预测；回归模型方法：大样本的内部预测；应用案例：应用案例：（1）CUMCM2004-A:奥运临时超市网点设计；奥运临时超市网点设计；（2）CUMCM2004-B:电力市场的输电阻塞管电力市场的输电阻塞管理；理；（3）CUMCM2005-A:长江水质的评价与预测；长江水质的评价与预测；（4）CUMCM2006-B:艾滋病疗法的评价与预艾滋病疗法的评价与预测；测；（5）CUMCM2008-B:高教学费标准探讨问题。高教学费标准探讨问题。1622020/12/2