首页 文章

贝叶斯与OLS

提问于
浏览
0

我在网上发现了这个问题 . 有人可以详细解释一下,为什么使用OLS更好?是不是因为样本数量不足?另外,为什么不使用所有1000个样本来估计先前的分布呢?

我们有1000个随机抽样的数据点 . 目标是尝试使用k个回归变量中的一个响应变量构建回归模型 . 哪个更好? 1.(贝叶斯回归)使用前500个样本来估计假设的先验分布的参数,然后使用最后500个样本来更新后验分布之前的后验估计,以用于最终回归模型 . 2.(OLS回归)使用包含所有1000个回归量变量的简单普通最小二乘回归模型

1 回答

  • 4

    “更好”始终是一个意见问题,它在很大程度上取决于背景 .

    Advantages to a frequentist OLS approach :更广泛,更快速,更容易被更广泛的受众访问(因此更难以解释) . 我的一位聪明的教授常说"You don't need to build an atom smasher when a flyswatter will do the trick."

    Advantages to an equivalent Bayesian approach :进一步模型开发更灵活,可以直接模拟派生/计算数量的后验(还有更多,但这些是我通过给定分析去贝叶斯的动机) . 注意“等价”这个词 - 你可以在贝叶斯框架中做一些你不能在频率论方法中做的事情 .

    嘿,这是R的探索,首先模拟数据,然后使用典型的OLS方法 .

    N <- 1000
    x <- 1:N
    epsilon <- rnorm(N, 0, 1)
    y <- x + epsilon
    
    summary(lm(y ~ x))
    ## 
    ## Call:
    ## lm(formula = y ~ x)
    ## 
    ## Residuals:
    ##     Min      1Q  Median      3Q     Max 
    ## -2.9053 -0.6723  0.0116  0.6937  3.7880 
    ## 
    ## Coefficients:
    ##              Estimate Std. Error  t value Pr(>|t|)    
    ## (Intercept) 0.0573955  0.0641910    0.894    0.371    
    ## x           0.9999997  0.0001111 9000.996   <2e-16 ***
    ## ---
    ## Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
    ## 
    ## Residual standard error: 1.014 on 998 degrees of freedom
    ## Multiple R-squared:      1,  Adjusted R-squared:      1 
    ## F-statistic: 8.102e+07 on 1 and 998 DF,  p-value: < 2.2e-16
    

    ...这里是一个等价的贝叶斯回归,在回归参数和所有1000个数据点上使用非信息先验 .

    library(R2jags)
    cat('model {
      for (i in 1:N){
        y[i] ~ dnorm(y.hat[i], tau)
        y.hat[i] <- a + b * x[i]
      }
      a ~ dnorm(0, .0001)
      b ~ dnorm(0, .0001)
      tau <- pow(sigma, -2)
      sigma ~ dunif(0, 100)
    }', file="test.jags")
    
    test.data <- list(x=x,y=y,N=1000)
    test.jags.out <- jags(model.file="test.jags", data=test.data, 
                      parameters.to.save=c("a","b","tau","sigma"), n.chains=3, n.iter=10000)
    test.jags.out$BUGSoutput$mean$a
    ## [1] 0.05842661
    test.jags.out$BUGSoutput$sd$a
    ## [1] 0.06606705
    test.jags.out$BUGSoutput$mean$b
    ## [1] 0.9999976
    test.jags.out$BUGSoutput$sd$b
    ## [1] 0.0001122533
    

    请注意,参数估计值和标准误差/标准偏差基本相同!

    现在这是另一个贝叶斯回归,使用前500个数据点来估计先验,然后使用前500个来估计后验 .

    test.data <- list(x=x[1:500],y=y[1:500],N=500)
    test.jags.out <- jags(model.file="test.jags", data=test.data, 
                      parameters.to.save=c("a","b","tau","sigma"), n.chains=3, n.iter=10000)
    
    cat('model {
      for (i in 1:N){
        y[i] ~ dnorm(y.hat[i], tau)
        y.hat[i] <- a + b * x[i]
      }
      a ~ dnorm(a_mn, a_prec)
      b ~ dnorm(b_mn, b_prec)
      a_prec <- pow(a_sd, -2)
      b_prec <- pow(b_sd, -2)
      tau <- pow(sigma, -2)
      sigma ~ dunif(0, 100)
    }', file="test.jags1")
    
    test.data1 <- list(x=x[501:1000],y=y[501:1000],N=500,
                       a_mn=test.jags.out$BUGSoutput$mean$a,a_sd=test.jags.out$BUGSoutput$sd$a,
                       b_mn=test.jags.out$BUGSoutput$mean$b,b_sd=test.jags.out$BUGSoutput$sd$b)
    test.jags.out1 <- jags(model.file="test.jags1", data=test.data1, 
                      parameters.to.save=c("a","b","tau","sigma"), n.chains=3, n.iter=10000)
    
    test.jags.out1$BUGSoutput$mean$a
    ## [1] 0.01491162
    test.jags.out1$BUGSoutput$sd$a
    ## [1] 0.08513474
    test.jags.out1$BUGSoutput$mean$b
    ## [1] 1.000054
    test.jags.out1$BUGSoutput$sd$b
    ## [1] 0.0001201778
    

    有趣的是,推论类似于OLS结果,但差不多 . 这让我怀疑用于训练先前的500个数据点在分析中没有像过去500那样多的重量,而之前的数据点实际上已经被淘汰了,尽管我不确定这一点 .

    无论如何,我想不出有什么理由不使用所有1000个数据点(和非信息先验),特别是因为我怀疑500 500使用前500和后500不同 .

    也许,所有这一切的答案是: I trust the OLS and 1000-point Bayesian results more than the 500+500, and OLS is simpler.

相关问题