回归分析

第三章 回归分析

3.1.1一元回归模型

一元一次回归主要使用公式计算，掌握公式以及写出代码即可

也可以使用系统提供代码解决，与多项式回归类似，不再赘述

3.1.2一元多项式回归

（1）多项式拟合

[p,S,mu]=polyfit(x,y,n)

p为回归系数向量，S是结构体，用于调用后续命令，mu是x均值与标准差数组（E（x），σ）

（2）多项式回归区间预测

[Y,Delta]=polyval(p,x0,S,alpha)

p、S是polyfit的输出，x0是自变量取值，alpha默认0.05，Y为预测值，Delta是Y置信区间半径

（3）多项式回归的GUI

polytool(x,y,n,alpha)

x,y为样本观测数据向量，n是多项式阶数（最高次数），alpha默认0.05

用于画出1.总体拟合图形2.(1-alpha)上、下置信区间的直线(红色)，左端数值框显示随自变量变化得到的预测数值以及1-alpha置信区间长度一半的数值

下面对一元二次案例进行详细分析

%例：为了分析X射线的杀菌作用，用200千伏的X射线来照射细菌，每次照射6分钟，平板计数法估计尚存活的细菌数，照射次数为t，照射后残留细菌为y
clear
t=1:15;
y=[352,211,197,160,142,106,104,60,56,38,36,32,21,19,15];
p=polyfit(t,y,2)     % 作二次多项式回归
% p=(1.99,-51.14,347.90)
% y=1.99t^2-51.14t+347.9
y1= polyval(p,t);    % 模型估计与作图
plot(t,y,'-*',t,y1,'-o');  
legend('原始数据','二次函数') 
xlabel('t(照射次数)'),  ylabel('y(残留细菌数)')
t0=16;
yc1= polyconf(p,t0)  % 预测t0=16时残留的细菌数
% yc1=39.04
polytool(t,y,2)		% 打开多项式回归的GUI

3.1.3一元非线性回归模型

（1）非线性拟合

fun=inline(’y’， ’参变量’， ’x’)

[beta,r,J]=nlinfit(x,y, ’model’，beta0)

x,y分别是n×m矩阵和n维列向量，对一元非线性回归，x为n维向量，model是事先定义的fun非线性函数，beta0是回归系数初值（解方程得到），beta是估计的最佳回归系数，r残差，J是Jacobin矩阵。

（2）非线性回归预测

ypred=nlpredci(fun,inputs,beta,r,J)

fun是拟合函数，inputs是要预测的自变量，beta、r、J是nlinfit的输出

（3）非线性回归置信区间

ci=nlparci(beta,r,J,alpha)

（4）非线性GUI

nlintool(x,y,fun,beta0)

下面针对钢包使用非线性双曲线回归模型拟合

拟合函数为y=x/(ax+b)

%例：炼钢厂出钢时所用盛钢水的钢包，由于钢水对耐火材料的侵蚀，容积不断增大，找出使用次数与增大容积之间的函数关系
clear
x=[2:16];
y=[6.42,8.2,9.58,9.5,9.7,10,9.93,9.99,10.49,10.59,10.6,10.8,10.6,10.9,10.76];
plot(x,y,'*')        %绘制散点图
%其次，确定回归系数的初值beta0。选择已知数据中的两点(2,6.42)和(16,10.76)代入
[a,b]=solve('6.42*(2*a+b)=2','10.76*(16*a+b)=16')
% 第三步，建立非线性双曲线回归模型
beta0=[0.084,0.1436];            %依第二步，初始参数值beta0(1)=a, beta0(2)=b
fun=inline('x./(b(1)*x+b(2))','b','x');  %建立内联函数 
[b,r,J]=nlinfit(x,y,fun,beta0);       %建立模型并求解最佳参数
b                              %输出最佳参数
hold on
y1=x./(b(1)*x+b(2));              %绘制拟合曲线
plot(x,y1,'-or')
legend('原始数据','拟合曲线')
xlabel('使用次数(x)'), ylabel('增大容积(y)')
%第四步，预测。在执行上面的程序中，继续输入命令:
ypred=nlpredci(fun,17,b,r,J)     %预测使用17次后的钢包容积10.96
%最后，确定回归模型参数的95%的置信区间
ci=nlparci(b,r,J)            %依据模型给出参数区间估计
nlintool(x,y,fun,beta0)		%打开GUI

对一元建模进行完整运用,代码换掉数据即可使用

% 例：在四川白鹅的生产性能研究中，得到如下一组关于雏鹅重x与70日龄重y的数据，建立回归方程，计算误差平方和及可决系数，x=85,95,115时预测值与区间
clear
x=[80 86 98 90 120 102 95 83 113 105 110 100]';   % 雏鹅重
y=[2350 2400 2720 2500 3150	2680 2630 2400 3080 2920 2960 2860]'; %70日龄重
plot(x,y,'*')                 %作散点图
xlabel('x(雏鹅重)')           %横坐标名
ylabel('y(70日龄重)')        %纵坐标名
%（2）建立直线回归方程。调用命令polyfit，从而求出参数的最小二乘估计。n= size(x,1)                      % 计算样本容量
[p,s]=polyfit(x,y,1);              % 调用命令polyfit计算回归参数
y1=polyval(p,x);                   % 计算回归模型的函数值
hold on
plot(x,y1)                         % 作回归方程的图形，结果如图3-15所示
p 
%（3）公式计算误差估计与可决系数。
TSS=sum((y-mean(y)).^2)              %计算总离差平方和
RSS=sum((y1-mean(y)).^2)             %计算回归平方和
ESS=sum((y-y1).^2)                   %计算残差平方和
R2=RSS/TSS                           %计算样本决定系数R2.
%（4）回归方程关系显著性的F检验
F=(n-2)*RSS/ESS               %计算的F统计量
F1=finv(0.95,1,n-2)           %查F统计量0.05的分位数
F2=finv(0.99,1,n-2)           %查F统计量0.01的分位数
%（5）回归关系显著性的t检验
T=p(2)/sqrt(ESS/(n-2))*sqrt(sum((x-mean(x)).^2))    %计算T统计量
T1=tinv(0.975,n-2)                             %t统计量0.05的分位数
T2=tinv(0.995,n-2)                             %t统计量0.01的分位数
%(6)预测
x1=[85,95,115]';                     % 输入自变量
yc=polyval(p,x1)                     % 计算预测值
[Y,Delta]=polyconf(p,x1,s);        
I1=[Y-Delta,Y+Delta]                % 预测置信区间
%在程序中加入：
polytool(x,y)                   % 交互功能
figure
bar(x,y-y1), legend('残差')        % 残差图 
h=lillietest(y-y1)                % 残差正态性检验