如何从 R 中的 ca.po 函数中提取 p 值?

How to extract p value from ca.po function in R?

我想获得两个 ca.po 模型的 p 值。有人可以告诉我怎么做吗?

在?

library(data.table)
library(urca)
dt_xy = as.data.table(timeSeries::LPP2005REC[, 2:3])
res = urca::ca.po(dt_xy, type = "Pz", demean = demean, lag = "short")
summary(res)

和结果。我在结果中标记了我需要的 p 值。

模型 1 p 值 = 0.9841

模型 2 p 值 = 0.1363

######################################## 
# Phillips and Ouliaris Unit Root Test # 
######################################## 

Test of type Pz 
detrending of series with constant and linear trend 

Response SPI :

Call:
lm(formula = SPI ~ zr + trd)

Residuals:
      Min        1Q    Median        3Q       Max 
-0.036601 -0.003494  0.000243  0.004139  0.024975 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)
(Intercept)  9.702e-04  7.954e-04   1.220    0.223
zrSPI       -1.185e-02  5.227e-02  -0.227    0.821
zrSII       -3.037e-02  1.374e-01  -0.221    0.825
trd         -6.961e-07  3.657e-06  -0.190    0.849

Residual standard error: 0.007675 on 372 degrees of freedom
Multiple R-squared:  0.0004236, Adjusted R-squared:  -0.007637 
F-statistic: 0.05255 on 3 and 372 DF,  p-value: 0.9841 **<--- I need this p.value**


Response SII :

Call:
lm(formula = SII ~ zr + trd)

Residuals:
       Min         1Q     Median         3Q        Max 
-0.0096931 -0.0018105 -0.0002734  0.0017166  0.0115427 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)  
(Intercept) -7.598e-05  3.012e-04  -0.252   0.8010  
zrSPI       -1.068e-02  1.979e-02  -0.540   0.5897  
zrSII       -9.574e-02  5.201e-02  -1.841   0.0664 .
trd          1.891e-06  1.385e-06   1.365   0.1730  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.002906 on 372 degrees of freedom
Multiple R-squared:  0.01476,   Adjusted R-squared:  0.006813 
F-statistic: 1.858 on 3 and 372 DF,  p-value: 0.1363 **<--- I need this p.value**



Value of test-statistic is: 857.4274 

Critical values of Pz are:
                  10pct    5pct     1pct
critical values 71.9586 81.3812 102.0167

您必须深入研究 res 对象才能查看其属性以及那里可用的内容。

attributes(reg)
...
#> 
#> $testreg
#> Response SPI :
#> 
#> Call:
#> lm(formula = SPI ~ zr + trd)
#> 
...

返回了一长串对象,但是我们可以看到在testreg下调用了lm的摘要,我们可以看到这是[=的属性之一14=]。我们也可以使用attr(res, "name")访问res的属性,所以让我们看一下testreg.

的组件
names(attributes(res))
#>  [1] "z"         "type"      "model"     "lag"       "cval"      "res"      
#>  [7] "teststat"  "testreg"   "test.name" "class"
names(attr(res, "testreg"))
#> [1] "Response SPI" "Response SII"

正如您在上面提到的,您正在寻找 2 个独立的 p 值,因为我们有两个独立的模型。让我们检索这些并看看它们是什么。

spi <- attr(res, "testreg")[["Response SPI"]]
sii <- attr(res, "testreg")[["Response SII"]]

class(spi)
#> [1] "summary.lm"

因此,它们每个都是一个 summary.lm 对象。有很多关于如何从 lmsummary.lm 对象中提取 p 值的文档,所以让我们使用 .

中描述的方法
get_pval <- function(summary_lm) {
  pf(
    summary_lm$fstatistic[1L],
    summary_lm$fstatistic[2L],
    summary_lm$fstatistic[3L],
    lower.tail = FALSE
  )
}

get_pval(spi)
#>     value 
#> 0.9840898
get_pval(sii)
#>     value 
#> 0.1363474

好了,这就是您感兴趣的两个 p 值!