时间戳

首页 > 解决方案 > Julia中回归系数的不同标准误差

r - Julia中回归系数的不同标准误差

问题描述

我正在搞乱标准错误，试图确保我理解它们以及我想要的，但看起来lm标准矩阵代数有不同的计算方法，或者我的计算可能不正确......

这是我跑的

using DataFrames, GLM, LinearAlgebra, Statistics, LaTeXStrings, PyPlot, Random, Distributions, ProgressBars

stderr = Array{Float64,2}(undef, 20,3)
manstderr = Array{Float64,2}(undef, 20,3)
sigstderr = Array{Float64,2}(undef, 20,3)

_stderr = Array{Float64,2}(undef, 20,2)
_manstderr = Array{Float64,2}(undef, 20,2)
_sigstderr = Array{Float64,2}(undef, 20,2)

σ = 10

for i in ProgressBar(1:Int(20))
    X = rand(Uniform(-100,100), Int(1e6), 2)
    y = 1X[:,1]+3X[:,2] + rand(Normal(0,σ),Int(1e6))
    data = DataFrame(x1=X[:,1],x2=X[:,2],y=y)
    
    # full
    ols = lm(@formula(y ~ x1 + x2), data)
    stderr[i,:] = stderror(ols)
    manstderr[i,:] = diag(var(y - predict(ols)) * inv(ols.mm.m'ols.mm.m))
    sigstderr[i,:] = diag(σ^2 * inv(ols.mm.m'ols.mm.m))
    
    # omit
    ols = lm(@formula(y ~ x1), data)
    _stderr[i,:] = stderror(ols)
    _manstderr[i,:] = diag(var(y - predict(ols)) * inv(ols.mm.m'ols.mm.m))
    _sigstderr[i,:] = diag(σ^2 * inv(ols.mm.m'ols.mm.m))
    
end
scatter(stderr[:,1],stderr[:,2], label="full")
scatter(_stderr[:,1],_stderr[:,2], label="omit")
_ = legend()

但是当我测量矩阵之间差异的范数时，它还不足以让我相信它们在做同样的事情。

println("full")
println(norm(stderr - manstderr),"\n",
        norm(stderr - sigstderr))
println()
println("omit")
println(norm(_stderr - _manstderr),"\n",
        norm(_stderr - _sigstderr))

full
0.044288628759881314
0.04428864660024369

omit
0.64138313684955
0.7755064744520012

如果它们以相同的方式计算，我希望规范应该更接近？我不认为这是由于舍入错误。

我认为应该manstderr如何lm找到标准错误，但我很惊讶它们并没有靠得更近。这就是为什么我也加入了sigstderr计算。

如果你很好奇，这个实验是为了看看不包括相关但不相关的回归量会如何影响标准误差。

我包含在R标签中是因为也许有些R人已经熟悉这里发生的事情。从我在其他帖子中读到的内容来看，似乎可能存在舍入错误，但我认为这不适用于我的情况。

标签： rjuliaregression

解决方案

推荐阅读