计量经济学简单线性回归OLS的Matlab程序
- 格式:pdf
- 大小:214.28 KB
- 文档页数:2
使用Matlab技术进行回归分析的基本步骤回归分析是统计学中一种用于研究变量间关系的方法,可以用来预测和解释变量之间的相关性。
在实际应用中,使用计算工具进行回归分析可以提高分析效率和准确性。
本文将介绍使用Matlab技术进行回归分析的基本步骤,并探讨其中的一些关键概念和技巧。
一、数据准备在进行回归分析之前,首先需要收集和整理相关的数据。
这些数据通常包括自变量和因变量。
自变量是用来解释或预测因变量的变量,而因变量是需要解释或预测的变量。
在Matlab中,可以将数据保存为数据矩阵,其中每一列代表一个变量。
二、模型建立在回归分析中,需要建立一个数学模型来描述自变量和因变量之间的关系。
最简单的线性回归模型可以表示为:Y = βX + ε,其中Y是因变量,X是自变量,β是回归系数,ε是误差项。
在Matlab中,可以使用regress函数来进行线性回归分析。
三、模型拟合模型拟合是回归分析的核心步骤,它的目标是找到最佳的回归系数,使得预测值与实际观测值之间的差异最小。
在Matlab中,可以使用OLS(Ordinary Least Squares)方法来进行最小二乘法回归分析。
该方法通过最小化残差平方和来估计回归系数。
四、模型诊断模型诊断是回归分析中非常重要的一步,它可以帮助我们评估模型的合理性和有效性。
在Matlab中,可以使用多种诊断方法来检验回归模型是否满足统计假设,例如残差分析、方差分析和假设检验等。
这些诊断方法可以帮助我们检测模型是否存在多重共线性、异方差性和离群值等问题。
五、模型应用完成模型拟合和诊断之后,我们可以使用回归模型进行一些实际应用。
例如,可以使用模型进行因变量的预测,或者对自变量的影响进行解释和分析。
在Matlab中,可以使用该模型计算新的观测值和预测值,并进行相关性分析。
六、模型改进回归分析并不是一次性的过程,我们经常需要不断改进模型以提高预测的准确性和解释的可靠性。
在Matlab中,可以使用变量选择算法和模型改进技术来优化回归模型。
Matlab 线性回归(拟合)对于多元线性回归模型:e x x y p p ++++=βββΛ110设变量12,,,p x x x y L 的n 组观测值为12(,,,)1,2,,i i ip i x x x y i n =L L .记 ⎪⎪⎪⎪⎪⎭⎫ ⎝⎛=np n n p p x x x x x x x x x x ΛΛΛΛΛΛΛΛ212222*********,⎪⎪⎪⎪⎪⎭⎫ ⎝⎛=n y y y y M 21,则⎪⎪⎪⎪⎪⎭⎫ ⎝⎛=p ββββM 10 的估计值为 y x x x b ')'(ˆ1-==β(11.2) 在Matlab 中,用regress 函数进行多元线性回归分析,应用方法如下:语法:b = regress(y, x)[b, bint, r, rint, stats] = regress(y, x)[b, bint, r, rint, stats] = regress(y, x, alpha)b = regress(y, x),得到的1+p 维列向量b 即为(11.2)式给出的回归系数β的估计值.[b, bint, r, rint, stats]=regress(y, x) 给出回归系数β的估计值b ,β的95%置信区间((1)2p +⨯向量)bint ,残差r 以及每个残差的95%置信区间(2⨯n 向量)rint ;向量stats 给出回归的R 2统计量和F 以及临界概率p 的值.如果i β的置信区间(bint 的第1i +行)不包含0,则在显著水平为α时拒绝0i β=的假设,认为变量i x 是显著的.[b, bint, r, rint, stats]=regress(y, x, alpha) 给出了bint 和rint 的100(1-alpha)%的置信区间.三次样条插值函数的MATLAB 程序matlab 的splinex = 0:10; y = sin(x); %插值点xx = 0:.25:10; %绘图点yy = spline(x,y,xx);plot(x,y,'o',xx,yy)非线性拟合非线性拟合可以用以下命令(同样适用于线形回归分析):1.beta = nlinfit(X,y,fun,beta0)X给定的自变量数据,Y给定的因变量数据,fun要拟合的函数模型(句柄函数或者内联函数形式),beta0函数模型中系数估计初值,beta返回拟合后的系数2.x = lsqcurvefit(fun,x0,xdata,ydata)fun要拟合的目标函数,x0目标函数中的系数估计初值,xdata自变量数据,ydata 函数值数据X拟合返回的系数(拟合结果)nlinfit格式:[beta,r,J]=nlinfit(x,y,’model’, beta0)Beta 估计出的回归系数r 残差J Jacobian矩阵x,y 输入数据x、y分别为n*m矩阵和n维列向量,对一元非线性回归,x为n维列向量。
线性回归代码实现(matlab)function J = computeCost(X, y, theta)%COMPUTECOST Compute cost for linear regression% J = COMPUTECOST(X, y, theta) computes the cost of using theta as the% parameter for linear regression to fit the data points in X and y% Initialize some useful valuesm = length(y); % number of training examples% You need to return the following variables correctlyJ = 0;% ====================== YOUR CODE HERE ======================% Instructions: Compute the cost of a particular choice of theta% You should set J to the cost.predictions = X * theta;sqrErrors = (predictions-y) .^ 2;J = 1/(2*m) * sum(sqrErrors);% =========================================================================end 1.1 详细解释转化成了向量(矩阵)形式,如果⽤其他的语⾔,⽤循环应该可以实现predictions = X * theta; % 这⾥的⼤X是矩阵sqrErrors = (predictions-y) .^ 2;function [theta, J_history] = gradientDescent(X, y, theta, alpha, num_iters)%GRADIENTDESCENT Performs gradient descent to learn theta% theta = GRADIENTDESENT(X, y, theta, alpha, num_iters) updates theta by% taking num_iters gradient steps with learning rate alpha% Initialize some useful valuesm = length(y); % number of training examplesJ_history = zeros(num_iters, 1);for iter = 1:num_iters% ====================== YOUR CODE HERE ====================== % Instructions: Perform a single gradient step on the parameter vector% theta.%% Hint: While debugging, it can be useful to print out the values% of the cost function (computeCost) and gradient here.%theta_temp = theta;for j = 1:size(X, 2)theta_temp(j) = theta(j)-alpha*(1/m)*(X*theta - y)' * X(:, j);endtheta = theta_temp;% ============================================================ % Save the cost J in every iterationJ_history(iter) = computeCost(X, y, theta);endend 2.1 解释J_history = zeros(num_iters, 1);theta_temp = theta; 把theta存起来。
用matlab做一元线性回归分析一元线性回归分析是在排除其他影响因素的假定其他影响因素确定的情况下,分析某一个因素(自变量)是如何影响另外一个事物(因变量)的过程,所进行的分析是比较理想化的。
用SPSS可以做一元线性回归分析,但是当回归的自变量比较多的时候,一个一个的输入会比较麻烦,增加了计算量,本文中描述了如何用matlab语言来实现一元线性回归分析。
在matlab中,regress命令是用来做回归的。
假如有96个SNP,作为自变量,有一个因变量,比如说HDL,LDL等等,将它们以列导入matlab。
值得注意的是:自变量前面必须有一列全为1的数据,看下面例子即可理解。
for i=1:96z=[ones(2334,1), x(:,i)];[b,bint,r,rint,stats]=regress(y,z);c(i,:)=stats;end在一元线性回归方程中,回归方程的显著性检验可以替代回归系数的显著性检验,并且F=T2百度中的一个例子:X=[1 1 4 6 8 11 14 17 21]'Y=[2.49 3.30 3.68 12.20 27.04 61.10 108.80 170.90 275.50]' X=[ones(9,1), X][b,bint,r,rint,stats]= regress(Y,X)输出向量b,bint为回归系数估计值和它们的置信区间,r,rint 为残差及其置信区间,stats是用于检验回归模型的统计量,有三个数值,第一个是R2,其中R是相关系数,第二个是F统计量值,第三个是与统计量F对应的概率P,当P<α时拒绝H0,回归模型成立regressMultiple linear regressionSyntaxb = regress(y,X)[b,bint] = regress(y,X)[b,bint,r] = regress(y,X)[b,bint,r,rint] = regress(y,X)[b,bint,r,rint,stats] = regress(y,X)[...] = regress(y,X,alpha)Descriptionb = regress(y,X) returns the least squares fit of y on X by solving the linear modelfo r β, where:y is an n-by-1 vector of observationsX is an n-by-p matrix of regressorsβ is a p-by-1 vector of parametersɛ is an n-by-1 vector of random disturbances[b,bint] = regress(y,X) returns a matrix bint of 95% confidence intervals for β.[b,bint,r] = regress(y,X) returns a vector, r of residuals.[b,bint,r,rint] = regress(y,X) returns a matrix rint of intervals that can be used to diagnose outliers. If rint(i,:) does not contain zero, then the ith residual is larger than would be expected, at the 5% significance level. This is evidence that the ith observation is an outlier.[b,bint,r,rint,stats] = regress(y,X) returns a vector stats that contains the R2 statistic, the F statistic and a p value for the full model, and an estimate of the error variance.。
利用 Matlab 作回归分析一元线性回归模型:2,(0,)y x N αβεεσ=++求得经验回归方程:ˆˆˆyx αβ=+ 统计量: 总偏差平方和:21()n i i SST y y ==-∑,其自由度为1T f n =-; 回归平方和:21ˆ()n i i SSR y y ==-∑,其自由度为1R f =; 残差平方和:21ˆ()n i i i SSE y y ==-∑,其自由度为2E f n =-;它们之间有关系:SST=SSR+SSE 。
一元回归分析的相关数学理论可以参见《概率论与数理统计教程》,下面仅以示例说明如何利用Matlab 作回归分析。
【例1】为了了解百货商店销售额x 与流通费率(反映商业活动的一个质量指标,指每元商品流转额所分摊的流通费用)y 之间的关系,收集了九个商店的有关数据,见下表1.试建立流通费率y 与销售额x 的回归方程。
表1 销售额与流通费率数据【分析】:首先绘制散点图以直观地选择拟合曲线,这项工作可结合相关专业领域的知识和经验进行,有时可能需要多种尝试。
选定目标函数后进行线性化变换,针对变换后的线性目标函数进行回归建模与评价,然后还原为非线性回归方程。
【Matlab数据处理】:【Step1】:绘制散点图以直观地选择拟合曲线x=[1.5 4.5 7.5 10.5 13.5 16.5 19.5 22.5 25.5];y=[7.0 4.8 3.6 3.1 2.7 2.5 2.4 2.3 2.2];plot(x,y,'-o')输出图形见图1。
510152025图1 销售额与流通费率数据散点图根据图1,初步判断应以幂函数曲线为拟合目标,即选择非线性回归模型,目标函数为:(0)b y ax b =< 其线性化变换公式为:ln ,ln v y u x == 线性函数为:ln v a bu =+【Step2】:线性化变换即线性回归建模(若选择为非线性模型)与模型评价% 线性化变换u=log(x)';v=log(y)';% 构造资本论观测值矩阵mu=[ones(length(u),1) u];alpha=0.05;% 线性回归计算[b,bint,r,rint,states]=regress(v,mu,alpha)输出结果:b =[ 2.1421; -0.4259]表示线性回归模型ln=+中:lna=2.1421,b=-0.4259;v a bu即拟合的线性回归模型为=-;y x2.14210.4259bint =[ 2.0614 2.2228; -0.4583 -0.3934]表示拟合系数lna和b的100(1-alpha)%的置信区间分别为:[2.0614 2.2228]和[-0.4583 -0.3934];r =[ -0.0235 0.0671 -0.0030 -0.0093 -0.0404 -0.0319 -0.0016 0.0168 0.0257]表示模型拟合残差向量;rint =[ -0.0700 0.02300.0202 0.1140-0.0873 0.0813-0.0939 0.0754-0.1154 0.0347-0.1095 0.0457-0.0837 0.0805-0.0621 0.0958-0.0493 0.1007]表示模型拟合残差的100(1-alpha)%的置信区间;states =[0.9928 963.5572 0.0000 0.0012] 表示包含20.9928SSR R SST==、 方差分析的F 统计量/963.5572//(2)R E SSR f SSR F SSE f SSE n ===-、 方差分析的显著性概率((1,2))0p P F n F =->≈; 模型方差的估计值2ˆ0.00122SSE n σ==-。
佛山科学技术学院上机报告课程名称数学应用软件上机项目 MATLAB统计工具箱中的回归分析命令专业班级一. 上机目的本节课我们认识了用MA TALB统计工具箱中的回归命令,主要有以下内容:regress命令即可用于多元回归分析也可用于一元线性回归,其格式如下:1.确定回归系数的命令是regress,用命令:b=regress(Y,X).2.求回归系数的点的估计和区间估计,并检验回归模型,用命令:[b,bint,r,rint,stats]=regress(Y,X,alpha)3.画出残差及其置性区间,用命令:rcoplot(r,rint)二元多项式回归:[p,S]=polyfit(x,y,2)二. 上机内容1.第十六章课后习题1;2.第十六章课后习题2;3.第十六章课后习题3。
三. 上机方法与步骤给出相应的问题分析及求解方法,并写出Matlab程序,并有上机程序显示截图。
第1题:要求一元线性回归方程及检验其显著性,用命令[b,bint,r,rint,stats]=regress(Y,X);求置信区间和预测值用命令rstool(x,y,'purequadratic')回归方程及检验其显著性:x=[20 25 30 35 40 45 50 55 60 65]';X=[ones(10,1) x];Y=[13.2 15.1 16.4 17.1 17.9 18.7 19.6 21.2 22.5 24.3]';[b,bint,r,rint,stats]=regress(Y,X);b,bint,stats残差分析,作残差图:rcoplot(r,rint)预测及作图:z=b(1)+b(2)*xplot(x,Y,'k+',x,z,'r')预测值及置信区间:x=[20 25 30 35 40 45 50 55 60 65]';y=[13.2 15.1 16.4 17.1 17.9 18.7 19.6 21.2 22.5 24.3]'; rstool(x,y,'purequadratic')第2题:要求二元多项式回归方程用命令[p,S]=polyfit(x,y,2)二元回归方程:用polyfit命令编程:x=0:2:20;y=[0.6 2.0 4.4 7.5 11.8 17.1 23.3 31.2 39.6 49.7 61.7]; [p,S]=polyfit(x,y,2)用regress命令编程:x=[0 2 4 6 8 10 12 14 16 18 20];Y=[0.6 2.0 4.4 7.5 11.8 17.1 23.3 31.2 39.6 49.7 61.7]'; X=[ones(11,1) x' (x.^2)'];[b,bint,r,rint,stats]=regress(Y,X);b,stats预测及作图:Y=polyconf(p,x,y)plot(x,y,'k+',x,Y,'r')第3题:要求lny a b x=+型回归方程,用命令[beta,r,J]=nlinfit(x',y','volum3',beta0);function y=volum3(beta,x)y=beta(1)+beta(2)*log(x);x=[2 3 4 5 7 9 12 14 17 21 28 56];y=[35 42 47 53 59 65 68 73 76 82 86 99];beta0=[1 1]';[beta,r,J]=nlinfit(x',y','volum3',beta0);beta四.上机结果学会了编写程序,运用上机语言求出问题结果,验证结果。
MATLAB线性回归二、一元线性回归2.1.命令 polyfit最小二乘多项式拟合[p,S]=polyfit(x,y,m)多项式y=a1xm+a2xm-1+…+amx+am+1其中x=(x1,x2,…,xm)x1…xm为(n*1)的矩阵;y为(n*1)的矩阵;p=(a1,a2,…,am+1)是多项式y=a1xm+a2xm-1+…+amx+am+1的系数;S是一个矩阵,用来估计预测误差。
2.2.命令 polyval多项式函数的预测值Y=polyval(p,x)求polyfit所得的回归多项式在x处的预测值Y;p是polyfit函数的返回值;x和polyfit函数的x值相同。
2.3.命令 polyconf 残差个案次序图[Y,DELTA]=polyconf(p,x,S,alpha)求polyfit所得的回归多项式在x处的预测值Y及预测值的显著性为1—alpha的置信区间DELTA;alpha缺省时为0.05。
p是polyfit函数的返回值;x和polyfit函数的x值相同;S和polyfit函数的S值相同。
2.4 命令 polytool(x,y,m)一元多项式回归命令2.5.命令regress多元线性回归(可用于一元线性回归)b=regress( Y, X )[b, bint,r,rint,stats]=regress(Y,X,alpha)b 回归系数bint 回归系数的区间估计r 残差rint 残差置信区间stats 用于检验回归模型的统计量,有三个数值:相关系数R2、F值、与F对应的概率p,相关系数R2越接近1,说明回归方程越显著;F 〉F1—α(k,n—k-1)时拒绝H0,F越大,说明回归方程越显著;与F对应的概率p 时拒绝H0,回归模型成立。
Y为n*1的矩阵;X为(ones(n,1),x1,…,xm)的矩阵;alpha显著性水平(缺省时为0。
05)。
三、多元线性回归3.1.命令 regress(见2.5)3.2.命令 rstool 多元二项式回归命令:rstool(x,y,’model’,alpha)x 为n*m矩阵y为 n维列向量model 由下列4个模型中选择1个(用字符串输入,缺省时为线性模型):linear(线性):purequadratic(纯二次):interaction(交叉):quadratic(完全二次):alpha 显著性水平(缺省时为0。
计量经济学简单线性回归OLS的Matlab程序
wxh1000
2011-09-21
先写OLS.m的M文件,用来代替regress函数;
(目前对regress函数不太了解,这里特别感谢潘晓炜同学的提醒)
-----------------------------------------------------------------------------------↓
function [beta_0 beta_1]=OLS(y,x)
%Ordinary Linear Regression
%其中x,y为样本构成的向量;
%回归方程为Simple regression: y=beta_0+x*beta_1+u;
%y_mean=mean(y);
%x_mean=mean(x);
%beta_1=((x-x_mean)*(y-y_mean)')/((x-x_mean)*(x-x_mean)');
%beta_0=y_mean-beta_1*x_mean;
%其中u为服从N(0,sigma^2)随机变量;
y_mean=mean(y);
x_mean=mean(x);
beta_1=((x-x_mean)*(y-y_mean)')/((x-x_mean)*(x-x_mean)');
beta_0=y_mean-beta_1*x_mean;
-----------------------------------------------------------------------------------↑
然后写OLS_test.m的M文件,用来进行模拟;
-----------------------------------------------------------------------------------↓
function [b_0 b_1]=OLS_test(beta_0,beta_1,n,a,b,sigma)
%已知beta_0,beta_1,由OLS回归得b_0,b_1.两者进行比较得到估计效果;
%y=beta_0+beta_1*x+u来得到;
%x为随机向量,u为服从N(0,sigma^2)随机变量;
%n为模拟数据量,比如1,10,100,1000等;
%x=a+b*rand(1,n);%产生(a,a+b)区间上的随机向量;
%mu= ;sigma= ;%随机矩阵服从均值为mu,方差为sigma的正态分布
%M= ;N= %M,N为产生[M,N]的随机矩阵
%x=mu+sqrt(sigma)*randn(M,N);%x为新生成的矩阵[M,N],服从均值为mu,方差为sigma的正态分布;
x=a+b*rand(1,n);%产生(a,a+b)区间上的随机向量;
%随机矩阵服从均值为0,方差为sigma的正态分布
u=sqrt(sigma)*randn(1,n);
y=beta_0+beta_1*x+u;
%用OLS函数进行回归即可:[beta_0 beta_1]=OLS(y,x);
[b_0 b_1]=OLS(y,x);
sprintf('已知参数\n\tbeta_0=%0.5g\n\tbeta_1=%0.5g\n模拟后,OLS估计值为\n\tbeta_0=%0.5g\n\tbeta_1=%0.5g',beta_0,beta_1,b_0,b_1) -----------------------------------------------------------------------------------↑
运行结果如下:
-----------------------------------------------------------------------------------↓
>> [b_0 b_1]=OLS_test(20,0.7,2,1,99,1)
ans =
已知参数
beta_0=20
beta_1=0.7
模拟后,OLS估计值为
beta_0=19.394
beta_1=0.72626
>> [b_0 b_1]=OLS_test(20,0.7,100,1,99,1)
ans =
已知参数
beta_0=20
beta_1=0.7
模拟后,OLS估计值为
beta_0=19.894
beta_1=0.70479
>> [b_0 b_1]=OLS_test(20,0.7,10000,1,99,1)
ans =
已知参数
beta_0=20
beta_1=0.7
模拟后,OLS估计值为
beta_0=19.988
beta_1=0.70018。