如何从 R 中的 rms 包的 ols() 模型中获取 p.values 系数

How to get p.values of coefficients from ols() model of rms package in R

我运行回归使用mtcars数据集和rms包的ols()函数如下:

> library(rms)
> mod_ols = ols(mpg~wt+qsec+am+hp+vs+drat+cyl, mtcars)                           
> 
> mod_ols

Linear Regression Model

ols(formula = mpg ~ wt + qsec + am + hp + vs + drat + cyl, data = mtcars)

                Model Likelihood     Discrimination    
                   Ratio Test           Indexes        
Obs       32    LR chi2     62.73    R2       0.859    
sigma 2.5703    d.f.            7    R2 adj   0.818    
d.f.      24    Pr(> chi2) 0.0000    g        6.519    

Residuals

    Min      1Q  Median      3Q     Max 
-3.3615 -1.7407 -0.2537  1.0893  4.5929 

          Coef    S.E.    t     Pr(>|t|)
Intercept 15.7714 16.5491  0.95 0.3501  
wt        -3.0826  1.0888 -2.83 0.0092  
qsec       0.7426  0.6632  1.12 0.2739  
am         2.5945  1.8802  1.38 0.1803  
hp        -0.0182  0.0162 -1.12 0.2728  
vs         0.2321  2.0291  0.11 0.9099  
drat       0.6387  1.5186  0.42 0.6778  
cyl        0.0309  0.9091  0.03 0.9732  

生成的模型结构如下:

> str(mod_ols)
List of 18
 $ coefficients     : Named num [1:8] 15.7714 -3.0826 0.7426 2.5945 -0.0182 ...
  ..- attr(*, "names")= chr [1:8] "Intercept" "wt" "qsec" "am" ...
 $ residuals        : Named num [1:32] -2.193 -1.823 -3.362 0.714 1.81 ...
  ..- attr(*, "names")= chr [1:32] "Mazda RX4" "Mazda RX4 Wag" "Datsun 710" "Hornet 4 Drive" ...
 $ effects          : Named num [1:32] -113.65 -29.12 -9.1 5.12 3.04 ...
  ..- attr(*, "names")= chr [1:32] "Intercept" "wt" "qsec" "am" ...
 $ rank             : int 8
 $ fitted.values    : Named num [1:32] 23.2 22.8 26.2 20.7 16.9 ...
  ..- attr(*, "names")= chr [1:32] "Mazda RX4" "Mazda RX4 Wag" "Datsun 710" "Hornet 4 Drive" ...
 $ assign           :List of 7
  ..$ wt  : int 2
  ..$ qsec: int 3
  ..$ am  : int 4
  ..$ hp  : int 5
  ..$ vs  : int 6
  ..$ drat: int 7
  ..$ cyl : int 8
 $ qr               :List of 5
  ..$ qr   : num [1:32, 1:8] -5.657 0.177 0.177 0.177 0.177 ...
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:32] "Mazda RX4" "Mazda RX4 Wag" "Datsun 710" "Hornet 4 Drive" ...
  .. .. ..$ : chr [1:8] "Intercept" "wt" "qsec" "am" ...
  .. ..- attr(*, "assign")= int [1:8] 0 1 2 3 4 5 6 7
  ..$ qraux: num [1:8] 1.18 1.05 1.08 1.18 1.03 ...
  ..$ pivot: int [1:8] 1 2 3 4 5 6 7 8
  ..$ tol  : num 0.0000001
  ..$ rank : int 8
  ..- attr(*, "class")= chr "qr"
 $ df.residual      : int 24
 $ var              : num [1:8, 1:8] 273.872 3.919 -9.225 -13.12 -0.054 ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:8] "Intercept" "wt" "qsec" "am" ...
  .. ..$ : chr [1:8] "Intercept" "wt" "qsec" "am" ...
 $ stats            : Named num [1:6] 32 62.733 7 0.859 6.519 ...
  ..- attr(*, "names")= chr [1:6] "n" "Model L.R." "d.f." "R2" ...
 $ linear.predictors: Named num [1:32] 23.2 22.8 26.2 20.7 16.9 ...
  ..- attr(*, "names")= chr [1:32] "Mazda RX4" "Mazda RX4 Wag" "Datsun 710" "Hornet 4 Drive" ...
 $ call             : language ols(formula = mpg ~ wt + qsec + am + hp + vs + drat + cyl, data = mtcars)
 $ terms            :Classes 'terms', 'formula' length 3 mpg ~ wt + qsec + am + hp + vs + drat + cyl
  .. ..- attr(*, "variables")= language list(mpg, wt, qsec, am, hp, vs, drat, cyl)
  .. ..- attr(*, "factors")= int [1:8, 1:7] 0 1 0 0 0 0 0 0 0 0 ...
  .. .. ..- attr(*, "dimnames")=List of 2
  .. .. .. ..$ : chr [1:8] "mpg" "wt" "qsec" "am" ...
  .. .. .. ..$ : chr [1:7] "wt" "qsec" "am" "hp" ...
  .. ..- attr(*, "term.labels")= chr [1:7] "wt" "qsec" "am" "hp" ...
  .. ..- attr(*, "order")= int [1:7] 1 1 1 1 1 1 1
  .. ..- attr(*, "intercept")= num 1
  .. ..- attr(*, "response")= int 1
  .. ..- attr(*, ".Environment")=<environment: R_GlobalEnv> 
  .. ..- attr(*, "predvars")= language list(mpg, wt, qsec, am, hp, vs, drat, cyl)
  .. ..- attr(*, "dataClasses")= Named chr [1:8] "numeric" "numeric" "numeric" "numeric" ...
  .. .. ..- attr(*, "names")= chr [1:8] "mpg" "wt" "qsec" "am" ...
  .. ..- attr(*, "formula")=Class 'formula' length 3 mpg ~ wt + qsec + am + hp + vs + drat + cyl
  .. .. .. ..- attr(*, ".Environment")=<environment: R_GlobalEnv> 
 $ Design           :List of 11
  ..$ name        : chr [1:7] "wt" "qsec" "am" "hp" ...
  ..$ label       : chr [1:7] "wt" "qsec" "am" "hp" ...
  ..$ units       : Named chr [1:7] "" "" "" "" ...
  .. ..- attr(*, "names")= chr [1:7] "wt" "qsec" "am" "hp" ...
  ..$ colnames    : chr [1:7] "wt" "qsec" "am" "hp" ...
  ..$ assume      : chr [1:7] "asis" "asis" "asis" "asis" ...
  ..$ assume.code : int [1:7] 1 1 1 1 1 1 1
  ..$ parms       : list()
  ..$ limits      : list()
  ..$ values      : list()
  ..$ nonlinear   :List of 7
  .. ..$ wt  : logi FALSE
  .. ..$ qsec: logi FALSE
  .. ..$ am  : logi FALSE
  .. ..$ hp  : logi FALSE
  .. ..$ vs  : logi FALSE
  .. ..$ drat: logi FALSE
  .. ..$ cyl : logi FALSE
  ..$ interactions: list()
 $ non.slopes       : num 1
 $ na.action        : NULL
 $ scale.pred       : chr "mpg"
 $ fail             : logi FALSE
 - attr(*, "class")= chr [1:3] "ols" "rms" "lm"
> 

但是,我无法在 mod_ols 结构中找到 p.values。如何从 rms 包的 ols() 函数生成的模型中获取系数的 P 值?感谢您的帮助。

使用summary.lm:

summary.lm(mod_ols)$coefficients[ , "Pr(>|t|)"]
#  Intercept          wt        qsec          am          hp          vs 
#0.350087341 0.009234825 0.273857810 0.180345204 0.272787068 0.909889126 
#       drat         cyl 
#0.677806300 0.973163043

以下每一行都可以解决问题:

summary.lm(mod_ols)$coefficients[ , 4] 
coef(summary.lm(mod_ols))[ , 4]
coef(summary.lm(mod_ols))[ , "Pr(>|t|)"]