精华内容
下载资源
问答
  • 多重共线性一般是在(1)时间序列数据和(2)横截面数据中会发生。 产生的影响 (1)OLS得到的回归参数估计值很不稳定 (2)回归系数的方差随共线性强度增加而增长 (3)系数的正负号得不到合理的解释 多重共线性的...

    1、多重共线性

    多重共线性一般是在(1)时间序列数据和(2)横截面数据中会发生。

    产生的影响

    (1)OLS得到的回归参数估计值很不稳定
    (2)回归系数的方差随共线性强度增加而增长
    (3)系数的正负号得不到合理的解释

    多重共线性的判定方法

    1、方差膨胀因子vif
    2、特征根判定法

    3、直观判定法:更加直观
    1.当增加或剔除一个自变量,或者改变一个观测值时,回归系数的估计值发生较大变化。
    2.从定性分析认为,一些重要的自变量在回归方程中没有通过显著性检验。
    3.有些自变量的回归系数所带正负号与定性分析结果违背。
    4.自变量的相关矩阵中,自变量间的相关系数较大。
    5.一些重要的自变量的回归系数的标准误差较大。

    消除多重共线性

    1、剔除一些不重要的解释变量
    2、增大样本容量
    3、回归系数的有偏估计
    统计学家还致力于改进古典的最小二乘法,提出以采用有偏估计为代价来提高估计量稳定性的方法,如:
    逐步回归
    岭回归法
    主成分回归法
    偏最小二乘法
    lasso回归
    适应性lasso回归

    2、糖尿病数据+逐步回归法

    先看一下数据:
    在这里插入图片描述
    原模型的方差膨胀因子:几乎都非常大
    在这里插入图片描述
    这是逐步回归法的:
    在这里插入图片描述
    去除了大部分变量,保留了部分变量,方差膨胀因子相对小很多。
    这是逐步回归法的残差与y的散点图:
    在这里插入图片描述
    以及模型:
    在这里插入图片描述
    大部分自变量通过了检验。
    代码:
    数据需要可以评论~

    library("lars")
    library("car")
    #file.choose()
    da<-read.csv("F:\\learning_kecheng\\huigui\\67章—多重共线性和lasso\\diabetes.csv")[,11:75]
    head(da)
    
    lm1<-lm(y~.,data = da)
    lm_step = step(lm1,trace = F)
    vif(lm1)
    vif = vif(lm_step)
    sort(vif(lm1),de=T)[1:5]
    sort(vif)
    summary(lm_step)
    
    xx = cor(da[,-1])
    kappa(xx,exact = TRUE)
    kappa(da[,-1],exact = TRUE)
    
    shapiro.test(lm_step$res)
    Anova(lm_step,type="III")	
    plot(da$y,lm_step$res)
    abline(h=0,lty=2)
    
    展开全文
  • 一、案例介绍 1、目的:利用上市公司当年的公开财务指标预测来年盈利情况最重要的投资人决策依据。 2、数据来源:随机抽取深市和沪市2002和2003年的500个上市公司样本预测来年的净资产收益率。 3、解释变量包括:...

    一、案例介绍

    1、目的:利用上市公司当年的公开财务指标预测来年盈利情况最重要的投资人决策依据。

    2、数据来源:随机抽取深市和沪市2002和2003年的500个上市公司样本预测来年的净资产收益率。

    3、解释变量包括:资产周转率、当年净资产收益率、债务资本比率、市盈率、应收账款/主营业务收入、主营业务利润、存货/资产总计(反映公司存货状况)、对数资产总计(反映公司规模)

    二、描述性分析

    1、各个标量的均值、最小值、中位数、最大数和标准差

    2、变量相关性分析:相关性矩阵

    3、当期净资产收益率和往期净资产收益率的散点图

    三、建立模型:

    1、多元线性回归模型:

    2、模型假设:

    (1)解释变量是非随机的,且各解释变量之间互不相关(多重共线性)

    (2)随机误差项具有零均值、同方差和不序列相关性

    (3)解释变量和随机项不相关

    (4)随机项满足正态分布

    总结即:随机项满足零均值、同方差、不序列相关的正态分布;解释变量和随机项不相关且解释变量之间互不相关

    3、参数估计:

    (1)最小二乘估计量:

    RSS=\sum (y_{i}-\hat{\beta_{0}}-\hat{\beta_{1}}x_{i1}- \hat{\beta_{2}}x_{i2}-...-\hat{\beta_{p}}x_{ip})^2

    (2)方差估计量:

    \hat{\sigma }^2=RSS/(n-p-1)

    (3)拟合优度:

    总平方和:SST=\sum (y_i-\bar{y})^2

    残差平方和:SSe=\sum (y_i-\bar{y})^2

    R-square:R^2=1-\frac{SSE}{SST}

    4、显著性检验:

    (1)F检验

    假设:H_0:\beta_i=0 vs H_1:\beta_i\neq 0

    检验统计量:F=\frac{(SST-SSE)/p)}{SSE/(n-p-1))}\sim F_{p,n-p-1}

    (2)t检验

    假设:H_0:\beta_i=0 vs H_1:\beta_i\neq 0

    检验统计量:T=\frac{\hat{\beta_i}}{\sqrt{\sigma ^2/n}}\sim t_{n-p-1}

    5、模型检验

    (1)异方差性

    (2)正态性检验:

    QQ图:残差的分位数和正态分布的分位数呈线性关系

    Shapiro-Wilk normality test

    Kolmogorov-Smirnov test

    (3)异常值检验:待补充

    Cook距离

    (4)多重共线性检验:

    见五介绍多重共线性

    四、变量选择与预测:

    只有三个变量显著性通过,但是无法排除其他变量是否有预测能力。从而我们通过AIC和BIC准则选择。原理:同时考虑到了模型复杂度和拟合效果。

    AIC=n(log(\frac{RSS}{n})+1+log(2\pi ))+2p

    BIC=n(log(\frac{RSS}{n})+1+log(2\pi))+logn*p

    五、多重共线性问题:

    1、变量相关性对模型造成的影响:

    (1)完全多重共线性会使OLS(普通最小二乘)系数矩阵方程 解不唯一(基本上不存在完全多重共线性,多是不完全多重共线性),不完全多重共线性会使OLS估计量的方差和标准误较大(因为),即使得估计精度很小和置信区间变宽。

    (2)多重共线性由于自变量之间的相关性,从而变量估计系数可能出现完全相反的符号或者难以置信的数值。

    (3)可能出现显著自变量回归系数不显著:因为标准误较大,从而t检验的t值较小,倾向于接受原假设。

    (4)R方值较高,但t值并不都是统计显著的。R²等于回归平方和在总平方和中所占的比率,即回归方程所能解释的因变量变异性的百分比。具体解释见补充资料1:回归拟合增加解释变量为什么增加拟合优度。方差膨胀因子越接近1,多重共线性越严重。这个时候R2越接近1。

    2、多重共线性的诊断方法:

    (1)R2较高但t值统计显著的不多。

    (2)解释变量两两高度相关。

    (3)方差膨胀因子

    3、方差膨胀因子:

    (1)考虑辅助回归:x_i=a+\sum_{j=1}^{n}b_jx_j+e

    (2)R_{i}^{2}是辅助回归的拟合优度

    (3)方差膨胀因子:VIF_i=\frac{1}{1-{R_{i}}^{2}}

    在一定程度上在多大程度上第i个变量所包含的信息被其他变量覆盖。一般认为小于10就没有多重共线性问题。

     

    展开全文
  • 运用随机模拟的方法构造数据比较广义逆和一般逆在求解最小二乘估计时的结果,并进行残差分析,比较两种方法在完全多重共线性和半完全多重共线性性中的优缺点,最后对进一步研究复多重共线性提出相应建议:在接下来的...
  • 案例背景介绍这是mei国50个州关于犯罪率的一组数据,包括人口、面积、收入、文盲率、高中毕业率、霜冻天数、犯罪率7个指标,现在我们想考察一下州犯罪率和其他因素间的关系。SPSS变量视图如下:研究目标是各州的...

    当只考察一个自变量对因变量的影响时,我们称之为简单一元线性回归,如果要多考察一些自变量,此时许多人习惯性将之称为多元线性回归,统计学上建议称之为多重线性回归,避免和多元统计方法冲突。

    案例背景介绍

    这是mei国50个州关于犯罪率的一组数据,包括人口、面积、收入、文盲率、高中毕业率、霜冻天数、犯罪率共7个指标,现在我们想考察一下州犯罪率和其他因素间的关系。SPSS变量视图如下:

    研究目标是各州的犯罪率(因变量),可能的因素(自变量)是人口、面积、收入、文盲率、高中毕业率、霜冻天数。因变量犯罪率连续数值变量,有多个自变量,从研究目标和数据类型来看,可选用多重线性回归分析。

    线性关系初步判断

    线性回归要求每个自变量和因变量之间存在线性关系,可以依靠相关分析和散点图来初步判断。

    犯罪率与文盲率、霜冻天数、高中毕业率、人口存在较为明显的线性关系,面积和其他变量普遍无关,越冷的地方文盲率越低、高中毕业率越高。

    有统计学意义的相关系数依次为:0.703(文盲率)、-0.539(霜冻天数)、-0.488(高中毕业率)、0.344(人口)。除因变量外其他因素两两间相关系数均在0.7以下,因素间没有强相关关系存在,初步提示共线性问题较弱。

    以上分析表明,并不是所有因素都有犯罪率存在明显线性关系,如果我们构建多重线性回归,这可能涉及到自变量筛选的问题,可优先选择逐步回归的方法。

    共线性问题是由于自变量间存在强相关关系造成的,它的存在对回归是有影响的,现在我们需要观察6个自变量间的共线性问题,最为常见的依据则是关注容忍度Tol和方差膨胀因子VIF。

    SPSS在线性回归中可以是输出这两个指标,来看一下具体情况:

    VIF是Tol的倒数,所以它们两个其实是一回事,我们只需要解读其一即可。一般认为如果某个自变量的容忍度Tol<0.1,则可能存在严重共线性问题。反过来就是VIF>10提示存在较为严重共线性问题。

    本例中所有自变量的Tol值大于0.2,提示没有特别严重的共线性问题,综合相关系数的表现我们说这组数据自变量间共线性问题并不严重可忽略。

    开始逐步线性回归

    线性回归还要求残差独立、残差正态性、残差方差齐次,这些内容我们可以在回归后做残差诊断,异常值影响也放在回归后进行检查判断。

    现在,我们开始逐步回归。

    在【统计】按钮对话框中,建议在默认选项上新增【共线性诊断】、残差【德宾沃森】、残差【个案诊断】(3倍标准差)。德宾沃森检验残差独立性,残差个案诊断排查离群点。

    除了考虑残差标准差检查离群点外,建议新增【库克距离】(它综合残差和杠杆值)来诊断强影响点。

    除了【残差直方图】【残差正态图】外,增加绘制一个以标准化预测值为横轴,标准化残差为纵轴的散点图,主要用来判断残差正态性、残差方差齐次基本条件。

    好了,其他参数默认设置。执行。

    回归分析结果解读

    逐步回归显示,6个自变量中的文盲率和人口数依次被纳入模型,其他自变量没有进入模型。前后两个模型我们依据调整后的R方,认为模型2更优秀,此时模型可解释因变量总变异的54.8%。不算高也不算低,还有待继续提升。

    德宾沃森统计量值为2.18,接近2,认为残差具有独立性,满足条件。

    此后我们就只读取模型2的结果。模型显著性检验,P<0.05,说明模型中人口数、文盲率至少有一个是有统计学意义的。模型有统计学意义。

    两个自变量对犯罪率影响均有统计学意义。共线性问题可忽略。先不着急写出方程式。

    标准化残差正态PP图,大多数点落在对角线直线上,可认为残差正态,满足条件。

    标准化残差散点图,各点分布相对均匀,没有发现可循规律变化 。可认为残差方差齐次。

    以上我们就残差独立性、残差正态性、残差方差齐次均作出诊断,认为均满足条件。

    接下来看看有没有异常值对模型拟合产生影响。

    首先看看标准化残差离群点。残差个案诊断表明Nevada州标化残差3.094>3,可能是一个离群点,可考虑处理。

    再看一下强影响点,最大COOK距离0.196<0.5,基本认为没有强影响点。

    回归分析结果优化

    我们现在尝试剔除Nevada州的个案数据,看看新的回归模型表现如何。

    调整后R方值=0.636,和上一个模型(0.548)相比,提升是明显的,即新模型拟合质量明显提升。

    现在写出多重线性回归方程式:

    Y=4.359*文盲率+0.000251*人口数+1.052

    本例中,文盲率的回归系数4.4,表示控制其他因素不变时,文盲率上升1%,犯罪率将会上升4.4%。总体来看,我们所得模型可解释各州犯罪率64%的方差(变异)。

    全文完

    图/文=数据小兵

    参考自:《R语言实战》第2版。

    好文推荐阅读

    本文配套案例数据下载

    加入博客配套知识星球,下载本案例数据文件,对照练习,有问题请在知识星球内讨论。

    数据小兵坚持写博客已经12年

    坚持写微信公号文章6年

    坚持更新SPSS视频课程2年

    坚持一对一答疑讨论2年

    欢迎加入SPSS视频课程

    竭诚服务

    展开全文
  • 案例3 线性回归之汽车贷款(代码)7 线性回归模型与诊断Step1、导入数据数据清洗Step2、相关性分析Step3、线性回归算法1、简单线性回归3、多元线性回归3.1 多元线性回归的变量筛选Step4、残差分析Step5、强影响点...

    7 线性回归模型与诊断

    数据说明:本数据是一份汽车贷款数据

    字段名 中文含义
    id id
    Acc 是否开卡(1=已开通)
    avg_exp 月均信用卡支出(元)
    avg_exp_ln 月均信用卡支出的自然对数
    gender 性别(男=1)
    Age 年龄
    Income 年收入(万元)
    Ownrent 是否自有住房(有=1;无=0)
    Selfempl 是否自谋职业(1=yes, 0=no)
    dist_home_val 所住小区房屋均价(万元)
    dist_avg_income 当地人均收入
    high_avg 高出当地平均收入
    edu_class 教育等级:小学及以下开通=0,中学=1,本科=2,研究生=3
    get_ipython().magic('matplotlib inline')
    
    import matplotlib.pyplot as plt
    import os
    import numpy as np
    import pandas as pd
    import statsmodels.api as sm
    from statsmodels.formula.api import ols
    

    Step1、导入数据和数据清洗

    (直接去除了缺失值)

    raw = pd.read_csv(r'.\data\creditcard_exp.csv', skipinitialspace=True)
    raw.head()
    
    id Acc avg_exp avg_exp_ln gender Age Income Ownrent Selfempl dist_home_val dist_avg_income age2 high_avg edu_class
    0 19 1 1217.03 7.104169 1 40 16.03515 1 1 99.93 15.932789 1600 0.102361 3
    1 5 1 1251.50 7.132098 1 32 15.84750 1 0 49.88 15.796316 1024 0.051184 2
    2 95 0 NaN NaN 1 36 8.40000 0 0 88.61 7.490000 1296 0.910000 1
    3 86 1 856.57 6.752936 1 41 11.47285 1 0 16.10 11.275632 1681 0.197218 3
    4 50 1 1321.83 7.186772 1 28 13.40915 1 0 100.39 13.346474 784 0.062676 2
    exp = raw[raw['avg_exp'].notnull()].copy().iloc[:, 2:].drop('age2',axis=1)
    
    exp_new = raw[raw['avg_exp'].isnull()].copy().iloc[:, 2:].drop('age2',axis=1)
    
    exp.describe(include='all')
    
    avg_exp avg_exp_ln gender Age Income Ownrent Selfempl dist_home_val dist_avg_income high_avg edu_class
    count 70.000000 70.000000 70.000000 70.000000 70.000000 70.000000 70.000000 70.000000 70.000000 70.000000 70.000000
    mean 983.655429 6.787787 0.285714 31.157143 7.424706 0.385714 0.028571 74.540857 8.005472 -0.580766 1.928571
    std 446.294237 0.476035 0.455016 7.206349 3.077986 0.490278 0.167802 36.949228 3.070744 0.432808 0.873464
    min 163.180000 5.094854 0.000000 20.000000 3.493900 0.000000 0.000000 13.130000 3.828842 -1.526850 0.000000
    25% 697.155000 6.547003 0.000000 26.000000 5.175662 0.000000 0.000000 49.302500 5.915553 -0.887981 1.000000
    50% 884.150000 6.784627 0.000000 30.000000 6.443525 0.000000 0.000000 65.660000 7.084184 -0.612068 2.000000
    75% 1229.585000 7.114415 1.000000 36.000000 8.494237 1.000000 0.000000 105.067500 9.123105 -0.302082 3.000000
    max 2430.030000 7.795659 1.000000 55.000000 16.900150 1.000000 1.000000 157.900000 18.427000 0.259337 3.000000

    Step2、相关性分析

    散点图

    exp.plot('Income', 'avg_exp', kind='scatter')
    plt.show()
    
    exp[['Income', 'avg_exp', 'Age', 'dist_home_val']].corr(method='pearson')
    
    Income avg_exp Age dist_home_val
    Income 1.000000 0.674011 0.369129 0.249153
    avg_exp 0.674011 1.000000 0.258478 0.319499
    Age 0.369129 0.258478 1.000000 0.109323
    dist_home_val 0.249153 0.319499 0.109323 1.000000

    Step3、线性回归算法

    1、简单线性回归

    lm_s = ols('avg_exp ~ Income', data=exp).fit()
    lm_s.summary()
    
    OLS Regression Results
    Dep. Variable: avg_exp R-squared: 0.454
    Model: OLS Adj. R-squared: 0.446
    Method: Least Squares F-statistic: 56.61
    Date: Thu, 06 Feb 2020 Prob (F-statistic): 1.60e-10
    Time: 10:14:45 Log-Likelihood: -504.69
    No. Observations: 70 AIC: 1013.
    Df Residuals: 68 BIC: 1018.
    Df Model: 1
    Covariance Type: nonrobust
    coef std err t P>|t| [0.025 0.975]
    Intercept 258.0495 104.290 2.474 0.016 49.942 466.157
    Income 97.7286 12.989 7.524 0.000 71.809 123.648
    Omnibus: 3.714 Durbin-Watson: 1.424
    Prob(Omnibus): 0.156 Jarque-Bera (JB): 3.507
    Skew: 0.485 Prob(JB): 0.173
    Kurtosis: 2.490 Cond. No. 21.4


    Warnings:
    [1] Standard Errors assume that the covariance matrix of the errors is correctly specified.
    # Predict-在原始数据集上得到预测值和残差
    pd.DataFrame([lm_s.predict(exp), lm_s.resid], index=['predict', 'resid']
                ).T.head()
    
    predict resid
    0 1825.141904 -608.111904
    1 1806.803136 -555.303136
    3 1379.274813 -522.704813
    4 1568.506658 -246.676658
    5 1238.281793 -422.251793
    # 在待预测数据集上得到预测值
    lm_s.predict(exp_new)[:5]
    
    2     1078.969552
    11     756.465245
    13     736.919530
    19     687.077955
    20     666.554953
    dtype: float64
    

    3、多元线性回归

    lm_m = ols('avg_exp ~ Age + Income + dist_home_val + dist_avg_income',
               data=exp).fit()
    lm_m.summary()
    
    OLS Regression Results
    Dep. Variable: avg_exp R-squared: 0.542
    Model: OLS Adj. R-squared: 0.513
    Method: Least Squares F-statistic: 19.20
    Date: Thu, 06 Feb 2020 Prob (F-statistic): 1.82e-10
    Time: 10:15:18 Log-Likelihood: -498.59
    No. Observations: 70 AIC: 1007.
    Df Residuals: 65 BIC: 1018.
    Df Model: 4
    Covariance Type: nonrobust
    coef std err t P>|t| [0.025 0.975]
    Intercept -32.0078 186.874 -0.171 0.865 -405.221 341.206
    Age 1.3723 5.605 0.245 0.807 -9.822 12.566
    Income -166.7204 87.607 -1.903 0.061 -341.684 8.243
    dist_home_val 1.5329 1.057 1.450 0.152 -0.578 3.644
    dist_avg_income 261.8827 87.807 2.982 0.004 86.521 437.245
    Omnibus: 5.234 Durbin-Watson: 1.582
    Prob(Omnibus): 0.073 Jarque-Bera (JB): 5.174
    Skew: 0.625 Prob(JB): 0.0752
    Kurtosis: 2.540 Cond. No. 459.


    Warnings:
    [1] Standard Errors assume that the covariance matrix of the errors is correctly specified.

    3.1 多元线性回归的变量筛选

    '''forward select'''
    def forward_select(data, response):
        remaining = set(data.columns)
        remaining.remove(response)
        selected = []
        current_score, best_new_score = float('inf'), float('inf')
        while remaining:
            aic_with_candidates=[]
            for candidate in remaining:
                formula = "{} ~ {}".format(
                    response,' + '.join(selected + [candidate]))
                aic = ols(formula=formula, data=data).fit().aic
                aic_with_candidates.append((aic, candidate))
            aic_with_candidates.sort(reverse=True)
            best_new_score, best_candidate=aic_with_candidates.pop()
            if current_score > best_new_score: 
                remaining.remove(best_candidate)
                selected.append(best_candidate)
                current_score = best_new_score
                print ('aic is {},continuing!'.format(current_score))
            else:        
                print ('forward selection over!')
                break
                
        formula = "{} ~ {} ".format(response,' + '.join(selected))
        print('final formula is {}'.format(formula))
        model = ols(formula=formula, data=data).fit()
        return(model)
    
    data_for_select = exp[['avg_exp', 'Income', 'Age', 'dist_home_val', 
                           'dist_avg_income']]
    lm_m = forward_select(data=data_for_select, response='avg_exp')
    print(lm_m.rsquared)
    
    aic is 1007.6801413968117,continuing!
    aic is 1005.4969816306302,continuing!
    aic is 1005.2487355956046,continuing!
    forward selection over!
    final formula is avg_exp ~ dist_avg_income + Income + dist_home_val 
    0.541151292841195
    

    Step4、残差分析

    ana1 = lm_s
    
    exp['Pred'] = ana1.predict(exp)
    exp['resid'] = ana1.resid
    exp.plot('Income', 'resid',kind='scatter')
    plt.show()
    
    # 遇到异方差情况,教科书上会介绍使用加权最小二乘法,但是实际上最常用的是对被解释变量取对数
    ana1 = ols('avg_exp ~ Income', data=exp).fit()
    ana1.summary()
    
    OLS Regression Results
    Dep. Variable: avg_exp R-squared: 0.454
    Model: OLS Adj. R-squared: 0.446
    Method: Least Squares F-statistic: 56.61
    Date: Thu, 06 Feb 2020 Prob (F-statistic): 1.60e-10
    Time: 10:15:54 Log-Likelihood: -504.69
    No. Observations: 70 AIC: 1013.
    Df Residuals: 68 BIC: 1018.
    Df Model: 1
    Covariance Type: nonrobust
    coef std err t P>|t| [0.025 0.975]
    Intercept 258.0495 104.290 2.474 0.016 49.942 466.157
    Income 97.7286 12.989 7.524 0.000 71.809 123.648
    Omnibus: 3.714 Durbin-Watson: 1.424
    Prob(Omnibus): 0.156 Jarque-Bera (JB): 3.507
    Skew: 0.485 Prob(JB): 0.173
    Kurtosis: 2.490 Cond. No. 21.4


    Warnings:
    [1] Standard Errors assume that the covariance matrix of the errors is correctly specified.
    ana2 = ols('avg_exp_ln ~ Income', exp).fit()
    exp['Pred'] = ana2.predict(exp)
    exp['resid'] = ana2.resid
    #exp.plot('Income', 'resid',kind='scatter')
    ana2.summary()
    
    OLS Regression Results
    Dep. Variable: avg_exp_ln R-squared: 0.403
    Model: OLS Adj. R-squared: 0.394
    Method: Least Squares F-statistic: 45.92
    Date: Thu, 06 Feb 2020 Prob (F-statistic): 3.58e-09
    Time: 10:15:59 Log-Likelihood: -28.804
    No. Observations: 70 AIC: 61.61
    Df Residuals: 68 BIC: 66.11
    Df Model: 1
    Covariance Type: nonrobust
    coef std err t P>|t| [0.025 0.975]
    Intercept 6.0587 0.116 52.077 0.000 5.827 6.291
    Income 0.0982 0.014 6.776 0.000 0.069 0.127
    Omnibus: 10.765 Durbin-Watson: 1.197
    Prob(Omnibus): 0.005 Jarque-Bera (JB): 12.708
    Skew: -0.688 Prob(JB): 0.00174
    Kurtosis: 4.569 Cond. No. 21.4


    Warnings:
    [1] Standard Errors assume that the covariance matrix of the errors is correctly specified.
    # 取对数会使模型更有解释意义(R-squre甚至减小,考虑因变量取对数)
    exp['Income_ln'] = np.log(exp['Income'])
    ana3 = ols('avg_exp_ln ~ Income_ln', exp).fit()
    exp['Pred'] = ana3.predict(exp)
    exp['resid'] = ana3.resid
    exp.plot('Income_ln', 'resid',kind='scatter')
    plt.show()
    ana3.summary()
    
    OLS Regression Results
    Dep. Variable: avg_exp_ln R-squared: 0.480
    Model: OLS Adj. R-squared: 0.473
    Method: Least Squares F-statistic: 62.87
    Date: Thu, 06 Feb 2020 Prob (F-statistic): 2.95e-11
    Time: 10:16:30 Log-Likelihood: -23.950
    No. Observations: 70 AIC: 51.90
    Df Residuals: 68 BIC: 56.40
    Df Model: 1
    Covariance Type: nonrobust
    coef std err t P>|t| [0.025 0.975]
    Intercept 5.0611 0.222 22.833 0.000 4.619 5.503
    Income_ln 0.8932 0.113 7.929 0.000 0.668 1.118
    Omnibus: 8.382 Durbin-Watson: 1.368
    Prob(Omnibus): 0.015 Jarque-Bera (JB): 8.074
    Skew: -0.668 Prob(JB): 0.0177
    Kurtosis: 3.992 Cond. No. 13.2


    Warnings:
    [1] Standard Errors assume that the covariance matrix of the errors is correctly specified.
    # 寻找最优的模型
    r_sq = {'exp~Income':ana1.rsquared, 'ln(exp)~Income':ana2.rsquared, 
            'ln(exp)~ln(Income)':ana3.rsquared}
    print(r_sq)
    
    {'exp~Income': 0.45429062315565294, 'ln(exp)~Income': 0.4030855555329651, 'ln(exp)~ln(Income)': 0.48039279938931057}
    

    Step5、强影响点分析

    # 法一:
    # Find outlier:
    exp['resid_t'] = (exp['resid'] - exp['resid'].mean()) / exp['resid'].std()
    exp[abs(exp['resid_t']) > 2]
    
    avg_exp avg_exp_ln gender Age Income Ownrent Selfempl dist_home_val dist_avg_income high_avg edu_class Pred resid Income_ln resid_t
    73 251.56 5.527682 0 29 5.1578 0 0 63.23 5.492947 -0.335147 0 6.526331 -0.998649 1.640510 -2.910292
    98 163.18 5.094854 0 22 3.8159 0 0 63.27 3.997789 -0.181889 0 6.257191 -1.162337 1.339177 -3.387317
    # Drop outlier
    exp2 = exp[abs(exp['resid_t']) <= 2].copy()
    ana4 = ols('avg_exp_ln ~ Income_ln', exp2).fit()
    exp2['Pred'] = ana4.predict(exp2)
    exp2['resid'] = ana4.resid
    exp2.plot('Income', 'resid', kind='scatter')
    plt.show()
    ana4.rsquared
    
    0.49397191385172456
    
    #法二:
    # statemodels包提供了更多强影响点判断指标
    from statsmodels.stats.outliers_influence import OLSInfluence
    OLSInfluence(ana3).summary_frame().head()
    
    dfb_Intercept dfb_Income_ln cooks_d standard_resid hat_diag dffits_internal student_resid dffits
    0 0.343729 -0.381393 0.085587 -1.319633 0.089498 -0.413732 -1.326996 -0.416040
    1 0.307196 -0.341294 0.069157 -1.201699 0.087409 -0.371907 -1.205702 -0.373146
    3 0.207619 -0.244956 0.044984 -1.440468 0.041557 -0.299947 -1.452165 -0.302382
    4 0.112301 -0.127713 0.010759 -0.575913 0.060926 -0.146693 -0.573062 -0.145967
    5 0.120572 -0.150924 0.022274 -1.221080 0.029011 -0.211064 -1.225579 -0.211842
    # ### 增加变量
    # 经过单变量线性回归的处理,我们基本对模型的性质有了一定的了解,接下来可以放入更多的连续型解释变量。在加入变量之前,要注意变量的函数形式转变。比如当地房屋均价、当地平均收入,其性质和个人收入一样,都需要取对数
    
    exp2['dist_home_val_ln'] = np.log(exp2['dist_home_val'])
    exp2['dist_avg_income_ln'] = np.log(exp2['dist_avg_income'])
    
    ana5 = ols('''avg_exp_ln ~ Age + Income_ln + 
               dist_home_val_ln + dist_avg_income_ln''', exp2).fit()
    exp2.plot('Income', 'resid', kind='scatter')
    plt.show()
    ana5.rsquared
    
    0.5529068646270383
    

    Step6、多重共线性分析(vif函数)

    Step regression is not always work.

    ana5.bse # The standard errors of the parameter estimates
    
    Intercept             0.317453
    Age                   0.005124
    Income_ln             0.568848
    dist_home_val_ln      0.058210
    dist_avg_income_ln    0.612197
    dtype: float64
    
    # The function "statsmodels.stats.outliers_influence.variance_inflation_factor" uses "OLS" to fit data, and it will generates a wrong rsquared. So define it ourselves!
    def vif(df, col_i):
        cols = list(df.columns)
        cols.remove(col_i)
        cols_noti = cols
        formula = col_i + '~' + '+'.join(cols_noti)
        r2 = ols(formula, df).fit().rsquared
        return 1. / (1. - r2)
    
    exog = exp2[['Income_ln', 'dist_home_val_ln',
                 'dist_avg_income_ln']]
    for i in exog.columns:
        print(i, '\t', vif(df=exog, col_i=i))
    
    
    Income_ln 	 36.653639058963186
    dist_home_val_ln 	 1.053596313570258
    dist_avg_income_ln 	 36.894876856102
    
    # Income_ln与dist_avg_income_ln具有共线性,使用“高出平均收入的比率”代替其中一个
    exp2['high_avg_ratio'] = exp2['high_avg'] / exp2['dist_avg_income']
    
    exog1 = exp2[['high_avg_ratio', 'dist_home_val_ln', 
                  'dist_avg_income_ln']]
    
    for i in exog1.columns:
        print(i, '\t', vif(df=exog1, col_i=i))
    
    high_avg_ratio 	 1.1230220802048871
    dist_home_val_ln 	 1.0527009887483532
    dist_avg_income_ln 	 1.1762825351755393
    
    var_select = exp2[['avg_exp_ln', 'high_avg_ratio', 
                       'dist_home_val_ln', 'dist_avg_income_ln']]
    ana7 = forward_select(data=var_select, response='avg_exp_ln')
    print(ana7.rsquared)
    
    aic is 23.816793700737392,continuing!
    aic is 20.830952279560805,continuing!
    forward selection over!
    final formula is avg_exp_ln ~ dist_avg_income_ln + dist_home_val_ln 
    0.552039773684598
    
    formula8 = '''
    avg_exp_ln ~ dist_avg_income_ln + dist_home_val_ln + 
    C(gender) + C(Ownrent) + C(Selfempl) + C(edu_class)
    '''
    ana8 = ols(formula8, exp2).fit()
    ana8.summary()
    
    OLS Regression Results
    Dep. Variable: avg_exp_ln R-squared: 0.873
    Model: OLS Adj. R-squared: 0.858
    Method: Least Squares F-statistic: 58.71
    Date: Thu, 06 Feb 2020 Prob (F-statistic): 1.75e-24
    Time: 11:17:40 Log-Likelihood: 35.337
    No. Observations: 68 AIC: -54.67
    Df Residuals: 60 BIC: -36.92
    Df Model: 7
    Covariance Type: nonrobust
    coef std err t P>|t| [0.025 0.975]
    Intercept 4.5520 0.212 21.471 0.000 4.128 4.976
    C(gender)[T.1] -0.4301 0.060 -7.200 0.000 -0.550 -0.311
    C(Ownrent)[T.1] 0.0184 0.045 0.413 0.681 -0.071 0.107
    C(Selfempl)[T.1] -0.3805 0.119 -3.210 0.002 -0.618 -0.143
    C(edu_class)[T.2] 0.2895 0.051 5.658 0.000 0.187 0.392
    C(edu_class)[T.3] 0.4686 0.060 7.867 0.000 0.349 0.588
    dist_avg_income_ln 0.9563 0.098 9.722 0.000 0.760 1.153
    dist_home_val_ln 0.0522 0.034 1.518 0.134 -0.017 0.121
    Omnibus: 3.788 Durbin-Watson: 2.129
    Prob(Omnibus): 0.150 Jarque-Bera (JB): 4.142
    Skew: 0.020 Prob(JB): 0.126
    Kurtosis: 4.208 Cond. No. 60.2


    Warnings:
    [1] Standard Errors assume that the covariance matrix of the errors is correctly specified.
    formula9 = '''
    avg_exp_ln ~ dist_avg_income_ln + dist_home_val_ln + 
    C(Selfempl) + C(gender):C(edu_class)
    '''
    ana9 = ols(formula9, exp2).fit()
    ana9.summary()
    
    OLS Regression Results
    Dep. Variable: avg_exp_ln R-squared: 0.914
    Model: OLS Adj. R-squared: 0.902
    Method: Least Squares F-statistic: 78.50
    Date: Thu, 06 Feb 2020 Prob (F-statistic): 1.42e-28
    Time: 11:17:48 Log-Likelihood: 48.743
    No. Observations: 68 AIC: -79.49
    Df Residuals: 59 BIC: -59.51
    Df Model: 8
    Covariance Type: nonrobust
    coef std err t P>|t| [0.025 0.975]
    Intercept 4.4098 0.178 24.839 0.000 4.055 4.765
    C(Selfempl)[T.1] -0.2945 0.101 -2.908 0.005 -0.497 -0.092
    C(edu_class)[T.2] 0.3164 0.045 7.012 0.000 0.226 0.407
    C(edu_class)[T.3] 0.5576 0.054 10.268 0.000 0.449 0.666
    C(gender)[T.1]:C(edu_class)[1] -0.0054 0.098 -0.055 0.956 -0.201 0.190
    C(gender)[T.1]:C(edu_class)[2] -0.4357 0.068 -6.374 0.000 -0.573 -0.299
    C(gender)[T.1]:C(edu_class)[3] -0.6001 0.065 -9.230 0.000 -0.730 -0.470
    dist_avg_income_ln 0.9893 0.078 12.700 0.000 0.833 1.145
    dist_home_val_ln 0.0654 0.029 2.278 0.026 0.008 0.123
    Omnibus: 5.023 Durbin-Watson: 1.722
    Prob(Omnibus): 0.081 Jarque-Bera (JB): 5.070
    Skew: -0.328 Prob(JB): 0.0793
    Kurtosis: 4.166 Cond. No. 61.2


    Warnings:
    [1] Standard Errors assume that the covariance matrix of the errors is correctly specified.

    Step7、正则算法

    1、岭回归

    # L1_wt参数为0则使用岭回归,为1使用lasso
    lmr = ols('avg_exp ~ Income + dist_home_val + dist_avg_income',
              data=exp).fit_regularized(alpha=1, L1_wt=0)
    
    lmr.summary()
    
    # ### LASSO算法
    lmr1 = ols('avg_exp ~ Age + Income + dist_home_val + dist_avg_income',
               data=exp).fit_regularized(alpha=1, L1_wt=1)
    lmr1.summary()
    

    2、使用scikit-learn进行正则化参数调优

    from sklearn.preprocessing import StandardScaler
    
    continuous_xcols = ['Age', 'Income', 'dist_home_val', 
                        'dist_avg_income']   #  抽取连续变量
    scaler = StandardScaler()  # 标准化
    X = scaler.fit_transform(exp[continuous_xcols])
    y = exp['avg_exp_ln']
    
    d:\Anaconda3\lib\site-packages\sklearn\preprocessing\data.py:645: DataConversionWarning: Data with input dtype int64, float64 were all converted to float64 by StandardScaler.
      return self.partial_fit(X, y)
    d:\Anaconda3\lib\site-packages\sklearn\base.py:464: DataConversionWarning: Data with input dtype int64, float64 were all converted to float64 by StandardScaler.
      return self.fit(X, **fit_params).transform(X)
    
    from sklearn.linear_model import RidgeCV
    
    alphas = np.logspace(-2, 3, 100, base=10)
    
    # Search the min MSE by CV
    rcv = RidgeCV(alphas=alphas, store_cv_values=True) 
    rcv.fit(X, y)
    
    RidgeCV(alphas=array([1.00000e-02, 1.12332e-02, ..., 8.90215e+02, 1.00000e+03]),
        cv=None, fit_intercept=True, gcv_mode=None, normalize=False,
        scoring=None, store_cv_values=True)
    
    print('The best alpha is {}'.format(rcv.alpha_))
    print('The r-square is {}'.format(rcv.score(X, y))) 
    # Default score is rsquared
    
    The best alpha is 0.2915053062825176
    The r-square is 0.47568267770194916
    
    X_new = scaler.transform(exp_new[continuous_xcols])
    np.exp(rcv.predict(X_new)[:5])
    
    d:\Anaconda3\lib\site-packages\ipykernel_launcher.py:1: DataConversionWarning: Data with input dtype int64, float64 were all converted to float64 by StandardScaler.
      """Entry point for launching an IPython kernel.
    
    array([759.67677561, 606.74024213, 661.20654568, 681.888929  ,
           641.06967182])
    
    cv_values = rcv.cv_values_
    n_fold, n_alphas = cv_values.shape
    
    cv_mean = cv_values.mean(axis=0)
    cv_std = cv_values.std(axis=0)
    ub = cv_mean + cv_std / np.sqrt(n_fold)
    lb = cv_mean - cv_std / np.sqrt(n_fold)
    
    plt.semilogx(alphas, cv_mean, label='mean_score')
    plt.fill_between(alphas, lb, ub, alpha=0.2)
    plt.xlabel("$\\alpha$")
    plt.ylabel("mean squared errors")
    plt.legend(loc="best")
    plt.show()
    
    # 手动选择正则化系数——根据业务判断
    
    # 岭迹图
    
    # In[42]:
    
    from sklearn.linear_model import Ridge
    
    ridge = Ridge()
    
    coefs = []
    for alpha in alphas:
        ridge.set_params(alpha=alpha)
        ridge.fit(X, y)
        coefs.append(ridge.coef_)
    
    
    # In[43]:
    
    ax = plt.gca()
    
    ax.plot(alphas, coefs)
    ax.set_xscale('log')
    plt.xlabel('alpha')
    plt.ylabel('weights')
    plt.title('Ridge coefficients as a function of the regularization')
    plt.axis('tight')
    plt.show()
    
    rcv.coef_
    
    array([ 0.03321449, -0.30956185,  0.05551208,  0.59067449])
    
    ridge.set_params(alpha=40)
    ridge.fit(X, y)
    ridge.coef_
    
    array([0.03293109, 0.09907747, 0.04976305, 0.12101456])
    
    ridge.score(X, y)
    
    0.4255673043353688
    
    np.exp(ridge.predict(X_new)[:5])
    
    array([934.79025945, 727.11042209, 703.88143602, 759.04342764,
           709.54172995])
    
    # lasso
    
    from sklearn.linear_model import LassoCV
    
    lasso_alphas = np.logspace(-3, 0, 100, base=10)
    lcv = LassoCV(alphas=lasso_alphas, cv=10) # Search the min MSE by CV
    lcv.fit(X, y)
    
    print('The best alpha is {}'.format(lcv.alpha_))
    print('The r-square is {}'.format(lcv.score(X, y))) 
    # Default score is rsquared
    
    The best alpha is 0.04037017258596556
    The r-square is 0.4426451069862233
    
    from sklearn.linear_model import Lasso
    
    lasso = Lasso()
    lasso_coefs = []
    for alpha in lasso_alphas:
        lasso.set_params(alpha=alpha)
        lasso.fit(X, y)
        lasso_coefs.append(lasso.coef_)
    ax = plt.gca()
    ax.plot(lasso_alphas, lasso_coefs)
    ax.set_xscale('log')
    plt.xlabel('alpha')
    plt.ylabel('weights')
    plt.title('Lasso coefficients as a function of the regularization')
    plt.axis('tight')
    plt.show()
    lcv.coef_
    
    array([0.        , 0.        , 0.02789489, 0.26549855])
    
    # 弹性网络
    from sklearn.linear_model import ElasticNetCV
    
    l1_ratio = [.1, .5, .7, .9, .95, .99, 1] #0取消
    
    encv = ElasticNetCV(l1_ratio=l1_ratio)
    encv.fit(X,y)
    
    print('The best l1_ratio is {}'.format(encv.l1_ratio_))
    print('The best alpha is {}'.format(encv.alpha_))
    
    d:\Anaconda3\lib\site-packages\sklearn\model_selection\_split.py:2053: FutureWarning: You should specify a value for 'cv' instead of relying on the default value. The default value will change from 3 to 5 in version 0.22.
      warnings.warn(CV_WARNING, FutureWarning)
    The best l1_ratio is 0.1
    The best alpha is 0.6293876843197391
    
    展开全文
  • 回归正则化方法(Lasso,Ridge和ElasticNet)在高维和数据集变量之间多重共线性情况下运行良好。   数学上,ElasticNet被定义为L1和L2正则化项的凸组合: 通过适当设置α,ElasticNet包含L1和L2正则化...
  • 线性回归的原理 应用于当Y值是连续变量的场合 公式 y=X*Beta ...X变量间不存在多重共线性;残差符合标准正态分布 了解数据 数据结构 Y变量定义 X变量类型 花费金额分布 分数据 INS:训练集(用来...
  • 关注一下~,更多商业数据分析案例等你来撩前言探索性数据分析、数据清洗与预处理和多元线性回归模型构建完毕后,为提升模型精度及其稳健性,还需进行许多操作。方差膨胀因子便是非常经典的一步,原理简单,实现优雅...
  • 关注一下~,更多商业数据分析案例等你来撩前言探索性数据分析、数据清洗与预处理和多元线性回归模型构建完毕后,为提升模型精度及其稳健性,还需进行许多操作。方差膨胀因子便是非常经典的一步,原理简单,实现优雅...
  • 经过探索分析发现该案例特征存在①异常值问题②多重共线性问题,通过尝试删除异常值、或不重要特征或有共线性特征来提高R方。本次练习共应用了线性回归、岭回归、Lasso回归算法,并最终确定用岭回归算法建立模型并...
  • 建立线性回归模型的完整步骤(附代码)

    千次阅读 多人点赞 2020-04-26 21:21:59
    目录零、案例简介一、数据诊断1 正态性检验1.1 检验方法1.1.1 直方图法1.1.2 PP图与QQ图1.1.3 Shapiro检验和K-S检验1.2 校正方法2 多重共线性检验2.1 检验方法2.1.1 方差膨胀因子VIF2.2 校正方法3 线性相关性检验3.1...
  • 回归入门及案例实战

    2020-02-25 17:18:34
    线性回归 普通线性回归 最小二乘法 线性:得名于f(x)=ax+b的图像的形象 很直观 就是...岭回归是加了二阶正则项的最小二乘,主要适用于过拟合严重或各变量之间存在多重共线性的时候,岭回归是有bias的,这里的bia...
  • 完整的R语言预测建模实例-从数据清理到建模预测

    万次阅读 多人点赞 2016-09-23 16:17:30
    概述本文使用Kaggle上的一个公开数据集,从数据导入,清理...(多重)共线性的挑战 预测因子的量纲差异 以上的几个主要挑战,对于熟悉机器学习的人来说,应该都是比较清楚的,这个案例中会涉及到五个挑战中的缺失值,量纲
  • 之前我们介绍了多元线性回归的原理, 又通过一个案例对多元线性回归模型进一步了解, 其中谈到自变量之间存在高度相关, 容易产生多重共线性问题, 对于多重共线性问题的解决方法有:删除自变量, 改变数据形式, 添加正则...
  • 统计学方法与数据分析(上下册)

    热门讨论 2013-12-29 11:32:47
    20章,分上、下两册,每册10章。各章均有大量习题。本书给出了大量的实际例子,这些例子涉及众多的学科和实际领域,但又不过于专门,容易理解。在大部分章节中都使用实例未引入主题,并把统计概念和这些非常实际的...
  • 《android 3d游戏开发技术详解与典型案例》分为两篇22章,第一篇以简单易懂的实例为依托,详细介绍了opengl es各方面的基础知识,第二篇则对7个真实案例的开发步骤进行了详细的介绍,逐步向读者讲解android 3d游戏...
  • 文章目录Lasso概念• 定义• Lasso处理多重共线性原理二、linear_model.Lasso 类案例:Lasso特征选取① 读取数据集② 划分训练集、测试集③ 对线性回归、岭回归、Lasso进行对比④ 学习曲线 Lasso概念 • 定义 LASSO...
  • R开发——预测模型

    2016-12-31 11:06:43
    交叉验证以及多个预测模型的比较全过程,注重在实际数据建模过程中的实际问题和挑战,主要包括以下五个方面的挑战:缺失值的挑战异常值的挑战不均衡分布的挑战(多重)共线性的挑战预测因子的量纲差异以上的几个主要...
  • 本书 20 章,分为四部分。第一部分为第 1 章内容,介绍了移动互联网与移动 Web 技术的相关知识,使读者对移动 Web 应用有一定的了解;第二部分为第 2 ~ 8 章,介绍了 HTML 5 各方面的知识点,重点介绍了绘图、...

空空如也

空空如也

1 2
收藏数 21
精华内容 8
热门标签
关键字:

多重共线性数据案例