利用MATLAB进行回归分析

1/19利用MATLAB进行回归分析一、实验目的：1．了解回归分析的基本原理，掌握MATLAB实现的方法；2.练习用回归分析解决实际问题。二、实验内容：题目1社会学家认为犯罪与收入低、失业及人口规模有关，对20个城市的犯罪率y（每10万人中犯罪的人数）与年收入低于5000美元家庭的百分比1x、失业率2x和人口总数3x（千人）进行了调查，结果如下表。序号y1x2x3x序号y1x2x3x111.216.56.25871114.518.16.07895213.420.56.46431226.923.17.4762340.726.39.36351315.719.15.8279345.316.55.36921436.224.78.6741524.819.27.312481518.118.66.5625612.716.55.96431628.924.98.3854720.920.26.419641714.917.96.7716835.721.37.615311825.822.48.692198.717.24.97131921.720.28.4595109.614.36.47492025.716.96.73353表12.26（1）若1x~3x中至多只许选择2个变量，最好的模型是什么？（2）包含3个自变量的模型比上面的模型好吗？确定最终模型。（3）对最终模型观察残差，有无异常点，若有，剔除后如何。理论分析与程序设计：为了能够有一个较直观的认识，我们可以先分别作出犯罪率y与年收入低于5000美元家庭的百分比1x、失业率2x和人口总数3x（千人）之间关系的散点图，根据大致分布粗略估计各因素造成的影响大小，再通过逐步回归法确定应该选择哪几个自变量作为模型。2/19编写程序如下：clc;clearall;y=[11.213.440.75.324.812.720.935.78.79.614.526.915.736.218.128.914.925.821.725.7];%犯罪率（人/十万人）x1=[16.520.526.316.519.216.520.221.317.214.318.123.119.124.718.624.917.922.420.216.9];%低收入家庭百分比x2=[6.26.49.35.37.35.96.47.64.96.46.07.45.88.66.58.36.78.68.46.7];%失业率x3=[587643635692124864319641531713749789576227937416258547169215953353];%总人口数（千人）figure(1),plot(x1,y,'*');figure(2),plot(x2,y,'*');figure(3),plot(x3,y,'*');X1=[x1',x2',x3'];stepwise(X1,y)运行结果与结论：3/19犯罪率与低收入散点图犯罪率与失业率散点图犯罪率与人口总数散点图4/19低收入与失业率作为自变量低收入与人口总数作为自变量5/19失业率与人口总数作为自变量在图中可以明显看出前两图的线性程度很好，而第三个图的线性程度较差，从这个角度来说我们应该以失业率和低收入为自变量建立模型。并且我们也可以从相关性角度来选取自变量，可以看出低收入与失业率作为自变量时的RMSE=4.64848；低收入与人口总数作为自变量时的RMSE=5.62245；失业率与人口总数作为自变量时的RMSE=5.04083。我们看到当低收入与失业率作为自变量时RMSE最小，因此如果选择两个变量作为自变量的会，它们是最适合的。并且可以得到三者的关系为：1234.07251.223934.39894yxx；对同时选取三个自变量的模型分析：如果我们将其三者同时选为自变量，我们发现RMSE=4.58978，比低收入与失业率二者作为自变量时稍微小了一点，不过我们也发现此时的X3系数为0.00076937，几乎为零，是可以忽略的，因此我们仍然选取两个自变量做最终的模型。关系函数仍为：1234.07251.223934.39894yxx6/19低收入、失业率与人口总数都作为自变量残差分析：对我们设定的最终模型运用残差分析，编写程序如下：clc;clearall;y=[11.213.440.75.324.812.720.935.78.79.614.526.915.736.218.128.914.925.821.725.7];%犯罪率（人/十万人）x1=[16.520.526.316.519.216.520.221.317.214.318.123.119.124.718.624.917.922.420.216.9];%低收入家庭百分比x2=[6.26.49.35.37.35.96.47.64.96.46.07.45.88.66.58.36.78.68.46.7];%失业率x3=[587643635692124864319641531713749789576227937416258547169215953353];%总人口数（千人）n=20;X2=[ones(n,1),x1',x2'];[b,bint,r,rint,s]=regress(y',X2);rcoplot(r,rint)运行结果如下：7/19我们应该剔除第18、20组数据，剔除后，运行源程序得到新的结果如下：8/19这时我们在重复本题开始时的做法，就可以得到最终的关系函数了。剔除不符数据后再次运行程序得到结果那么最终的函数关系便为：1235.70951.602283.39259yxx简要分析：从最终得到的结果上来看，失业率与低收入都将导致犯罪的上升。通过本道例题让我们学会运用逐步回归命令stepwise来分析多自变量情况下的最优模型问题，得到最优模型后，我们再运用残差法找到不符的数据，将其剔除，这样我们就会得到一个比较科学准确的关系式，这个思路对我们分析回归问题很有效。题目2一家洗衣粉制造公司新产品实验时，关心洗衣粉泡沫的高度y与搅拌程度X1和洗衣粉用量X2之间的关系，其中搅拌程度从弱到强分为3个水平。实验得到的数据如下表：X1X2Y1628.11732.31834.81938.29/1911043.5X1X2Y2665.32767.72869.42972.221076.9X1X2Y3682.23785.33888.13990.731093.6表12.30（1）将搅拌程度X1作为普通变量，建立y与X1和X2的回归模型，从残差图上发现问题。（2）将搅拌程度X1视为没有定量关系的3个水平，用0-1变量表示，建立回归模型，与（1）比较，从残差图上还能发现什么问题。（3）加入搅拌程度与洗衣粉用量的交互项，看看模型有无改进。理论分析与程序设计：仿照题目1中的程序，我们对搅拌程度（当成普通变量）与洗衣粉用量建立回归模型，并且进行残差分析。编写程序如下：clc;clearall;y=[28.132.334.838.243.565.367.769.472.276.982.285.388.190.793.6];%洗衣粉泡沫高度x1=[111112222233333];%搅拌程度x2=[678910678910678910];%洗衣粉用量figure(1),plot(x1,y,'*');figure(2),plot(x2,y,'*');X1=[x1',x2'];stepwise(X1,y)10/19运行结果如下：搅拌程度与泡沫高度关系洗衣粉用量与泡沫高度的关系11/19我们还可以得到选取不同的变量时RMSE的大小，为了让RMSE最小，我们根据“Nextstep”的提示，最后得到如下结果：从上图中可以看出当含有x1、x2两项时RMSE最小，因此该模型建立为：1212.7426.33.08667yxx编写程序进行残差分析：clc;clearall;y=[28.132.334.838.243.565.367.769.472.276.982.285.388.190.793.6];%洗衣粉泡沫高度x1=[111112222233333];%搅拌程度x2=[678910678910678910];%洗衣粉用量figure(1),plot(x1,y,'*');figure(2),plot(x2,y,'*');X1=[x1',x2'];stepwise(X1,y)n=15;X2=[ones(n,1),x1',x2'];[b,bint,r,rint,s]=regress(y',X2);rcoplot(r,rint)我们得到如下结果：12/19说明中等搅拌程度的残差与其他不同。将搅拌程度用0-1表示，泡沫高度关系分析：我们不妨设x0与x1共同表示搅拌程度，当（x0，x1）=（0,0）时，代表搅拌程度为1；当（x0，x1）=（0,1）时，代表搅拌程度为2；当（x0，x1）=（1,0）时，代表搅拌程度为3，我们就可以编写如下程序：clc;clearall;y=[28.132.334.838.243.565.367.769.472.276.982.285.388.190.793.6];%洗衣粉泡沫高度x0=[000000000011111];x1=[000001111100000];%搅拌程度x2=[678910678910678910];%洗衣粉用量figure(1),plot(x0,y,'*');figure(2),plot(x1,y,'*');figure(3),plot(x2,y,'*');X1=[x0',x1',x2'];stepwise(X1,y)n=15;X2=[ones(n,1),x0',x1',x2'];13/19[b,bint,r,rint,s]=regress(y',X2);rcoplot(r,rint)运行结果与分析：Y与x1关系散点图14/19Y与x2关系散点图Y与x3关系散点图15/19将搅拌程度改为0-1后的运行结果可以看出当RMSE最小时含有x0、x1、x2项，此时的函数关系式为：01210.686752.634.923.08667yxxx接下来我们看一下残差分析结果：16/19我们发现第5组数据应该剔除，剔除后再次运行程序，得到如下结果：而这时回归结果显示：这时的函数关系为：01211.6653.18435.5042.892yxxx17/19将交互项考虑进来分析：改写程序如下：clc;clearall;y=[28.132.334.838.243.565.367.769.472.276.982.285.388.190.793.6];%洗衣粉泡沫高度x1=[111112222233333];%搅拌程度x2=[678910678910678910];%洗衣粉用量figure(1),plot(x1,y,'*');figure(2),plot(x2,y,'*');X1=[x1',x2',(x1.*x2)'];stepwise(X1,y)n=15;X2=[ones(n,1),x1',x2',(x1.*x2)'];[b,bint,r,rint,s]=regress(y',X2);rcoplot(r,rint)clc;clearall;y=[28.132.334.838.243.565.367.769.472.276.982.285.388.190.793.6];%洗衣粉泡沫高度x0=[000000000011111];x1=[000001111100000];%搅拌程度x2=[678910678910678910];%洗衣粉用量figure(1),plot(x0,y,'*');figure(2),plot(x1,y,'*');figure(3),plot(x2,y,'*');X1=[x0',x1',x2',(x1.*x2)',(x1.*x0)',(x0.*x2)'];stepwise(X1,y)n=15;X2=[ones(n,1),x0',x1',x2',(x1.*x2)',(x1.*x0)',(x0.*x2)'