当我在R中拟合多项式回归模型时,如果我使用函数poly(),然后尝试用vif()求出方差通货膨胀因子,则得到以下误差:
y = c(0.22200,0.39500,0.42200,0.43700,0.42800,0.46700,0.44400,0.37800,0.49400,
0.45600,0.45200,0.11200,0.43200,0.10100,0.23200,0.30600,0.09230,0.11600,
0.07640,0.43900,0.09440,0.11700,0.07260,0.04120,0.25100,0.00002)
x1 = c(7.3,8.7,8.8,8.1,9.0,8.7,9.3,7.6,10.0,8.4,9.3,7.7,9.8,7.3,8.5,9.5,7.4,7.8,
7.7,10.3,7.8,7.1,7.7,7.4,7.3,7.6)
x2 = c(0.0,0.0,0.7,4.0,0.5,1.5,2.1,5.1,0.0,3.7,3.6,2.8,4.2,2.5,2.0,2.5,2.8,2.8,
3.0,1.7,3.3,3.9,4.3,6.0,2.0,7.8)
x3 = c(0.0,0.3,1.0,0.2,1.0,2.8,1.0,3.4,0.3,4.1,2.0,7.1,2.0,6.8,6.6,5.0,7.8,7.7,
8.0,4.2,8.5,6.6,9.5,10.9,5.2,20.7)
m = lm(y~poly(x1, x2, x3, degree=2, raw=TRUE))
summary(m)现在打电话给vif()
> vif(m)
Error in vif.default(m) : model contains fewer than 2 terms该模型具有9项和一个截距。
> m$rank
[1] 10在我看来,vif()函数不适用于poly()。这是正确的吗?有办法解决这个问题吗?还是我需要使用基本原则?
我可以计算出如下的方差通货膨胀因素:
X = poly(x1, x2, x3, degree=2, raw=TRUE)
C = solve(t(X)%*%X)
vifs = 1/diag(C)发布于 2022-11-10 22:16:06
我认为将poly与raw=TRUE结合使用没有多大意义,尤其是在这种情况下,您需要多个术语,而poly提供了一个非常不透明的结果标签。目前还不清楚使用的是哪个版本的vif。我选择了使用I()来创建“纯”二级术语的实验(.“同质”在这里是正确的吗?)和R公式接口(使用(...)^2)为其余,没有困难了解的结果,与我得到的与poly;
> m = lm(y~(x1+x2+x3)^2 + I(x1^2) + I(x2^2) + I(x3^2))
> summary(m)
Call:
lm(formula = y ~ (x1 + x2 + x3)^2 + I(x1^2) + I(x2^2) + I(x3^2))
Residuals:
Min 1Q Median 3Q Max
-0.063213 -0.037282 -0.001113 0.016738 0.122539
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -1.769364 1.286976 -1.375 0.1881
x1 0.420799 0.294173 1.430 0.1718
x2 0.222453 0.130742 1.701 0.1082
x3 -0.127995 0.070245 -1.822 0.0872 .
I(x1^2) -0.019325 0.016797 -1.150 0.2668
I(x2^2) -0.007449 0.012048 -0.618 0.5451
I(x3^2) 0.000824 0.001441 0.572 0.5754
x1:x2 -0.019876 0.012037 -1.651 0.1182
x1:x3 0.009151 0.007621 1.201 0.2473
x2:x3 0.002576 0.007039 0.366 0.7192
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.06092 on 16 degrees of freedom
Multiple R-squared: 0.9169, Adjusted R-squared: 0.8702
F-statistic: 19.63 on 9 and 16 DF, p-value: 5.051e-07
> rms::vif(m)
x1 x2 x3 I(x1^2) I(x2^2) I(x3^2) x1:x2 x1:x3 x2:x3
521.01297 401.58833 688.02220 501.50614 173.60055 99.67708 204.43081 456.00750 349.97018 为了找出vif的哪个版本“不起作用”,我尝试了rms::vif和HH::vif,并且没有一个抛出您遇到的错误,所以我不知道为什么会出现错误:
> m = lm(y~poly(x1, x2, x3, degree=2, raw=TRUE))
> HH::vif(m)
poly(x1, x2, x3, degree = 2, raw = TRUE)1.0.0 poly(x1, x2, x3, degree = 2, raw = TRUE)2.0.0
521.01297 501.50614
poly(x1, x2, x3, degree = 2, raw = TRUE)0.1.0 poly(x1, x2, x3, degree = 2, raw = TRUE)1.1.0
401.58833 204.43081
poly(x1, x2, x3, degree = 2, raw = TRUE)0.2.0 poly(x1, x2, x3, degree = 2, raw = TRUE)0.0.1
173.60055 688.02220
poly(x1, x2, x3, degree = 2, raw = TRUE)1.0.1 poly(x1, x2, x3, degree = 2, raw = TRUE)0.1.1
456.00750 349.97018
poly(x1, x2, x3, degree = 2, raw = TRUE)0.0.2
99.67708
> rms::vif(m)
poly(x1, x2, x3, degree = 2, raw = TRUE)1.0.0 poly(x1, x2, x3, degree = 2, raw = TRUE)2.0.0
521.01297 501.50614
poly(x1, x2, x3, degree = 2, raw = TRUE)0.1.0 poly(x1, x2, x3, degree = 2, raw = TRUE)1.1.0
401.58833 204.43081
poly(x1, x2, x3, degree = 2, raw = TRUE)0.2.0 poly(x1, x2, x3, degree = 2, raw = TRUE)0.0.1
173.60055 688.02220
poly(x1, x2, x3, degree = 2, raw = TRUE)1.0.1 poly(x1, x2, x3, degree = 2, raw = TRUE)0.1.1
456.00750 349.97018
poly(x1, x2, x3, degree = 2, raw = TRUE)0.0.2
99.67708 也许是car送的
> car::vif(m)
Error in vif.default(m) : model contains fewer than 2 terms如果是这样的话,正确的处理方法似乎是将问题提交给约翰·福克斯。他有一个GitHub页面来解决问题吗?通过查看输出packageDescription("car"),您应该能够找到这个问题的答案。不是的。你需要给他发封电子邮件。以下是如何到达目的地的方法:
maintainer("car")增编:
m = lm(y~(x1+x2+x3)^2 + I(x1^2) + I(x2^2) + I(x3^2))
> car::vif(m)
there are higher-order terms (interactions) in this model
consider setting type = 'predictor'; see ?vif
x1 x2 x3 I(x1^2) I(x2^2) I(x3^2) x1:x2 x1:x3 x2:x3
521.01297 401.58833 688.02220 501.50614 173.60055 99.67708 204.43081 456.00750 349.97018 所以也许在帮助页面中有一些有用的信息?也许值得一看这条信息到底是什么意思。
https://stackoverflow.com/questions/74394666
复制相似问题