Extraer pvalue de glm

Question

Extraer pvalue de glm

Estoy ejecutando muchas regresiones y solo estoy interesado en el efecto sobre el coeficiente y el valor p de una variable en particular. Por lo tanto, en mi script, me gustaría poder extraer el valor p del resumen glm (obtener el coeficiente en sí es fácil). La única forma que conozco de ver el valor p es usando summary (myReg). ¿Hay alguna otra manera?

Ej:

fit <- glm(y ~ x1 + x2, myData)
x1Coeff <- fit$coefficients[2] # only returns coefficient, of course
x1pValue <- ???

He intentado tratar fit$coefficients como una matriz, pero todavía no puedo simplemente extraer el valor p.

Es es posible hacer esto?

Gracias!

27

r glm p-value

Author: Clark Henry, 2014-05-24

Source

4 answers

He utilizado esta técnica en el pasado para extraer datos predictores de summary o de un objeto modelo ajustado:

coef(summary(m))[grepl("var_i_want$",row.names(coef(summary(m)))), 4]

Que me permite editar fácilmente en qué variable quiero obtener datos.

O como se señaló ser @ Ben, use match o %in%, algo más limpio que grepl:

coef(summary(m))[row.names(coef(summary(m))) %in% "var_i_want" , 4]

4

Author: James King,
Warning: date(): Invalid date.timezone value 'Europe/Kyiv', we selected the timezone 'UTC' for now. in /var/www/agent_stack/data/www/ajaxhispano.com/template/agent.layouts/content.php on line 61
2014-05-23 22:23:41

En lugar del número puede poner directamente el nombre

coef(summary(fit))[,'Pr(>|z|)']

Los otros disponibles en el resumen del coeficiente:

Estimate Std. Error z value Pr(>|z|)

1

Author: R. Prost,
Warning: date(): Invalid date.timezone value 'Europe/Kyiv', we selected the timezone 'UTC' for now. in /var/www/agent_stack/data/www/ajaxhispano.com/template/agent.layouts/content.php on line 61
2018-02-07 15:06:52

Bueno, esta sería otra manera, sin embargo, no la forma más eficiente de realizarla:

a = coeftable(model).cols[4]
pVals = [ a[i].v for i in 1:length(a) ]

Esto asegura que los valores extraídos del glm no estén en StatsBase. En él, puede jugar con pVals según su requisito. Espero que ayude, Ebby

1

Author: Ebby,
Warning: date(): Invalid date.timezone value 'Europe/Kyiv', we selected the timezone 'UTC' for now. in /var/www/agent_stack/data/www/ajaxhispano.com/template/agent.layouts/content.php on line 61
2018-05-20 11:48:59

score 42 · Accepted Answer

Quieres

coef(summary(fit))[,4]

Que extrae el vector de columna de p valores de la salida tabular mostrada por summary(fit). Los p-valores no se calculan realmente hasta que se ejecuta summary() en el ajuste del modelo.

Por cierto, use funciones extractor en lugar de profundizar en objetos si puede:

fit$coefficients[2]

Debe ser

coef(fit)[2]

Si no hay funciones extractor, str() es tu amigo. Le permite mirar la estructura de cualquier objeto, lo que le permite vea lo que contiene el objeto y cómo extraerlo:

summ <- summary(fit)

> str(summ, max = 1)
List of 17
 $ call          : language glm(formula = counts ~ outcome + treatment, family = poisson())
 $ terms         :Classes 'terms', 'formula' length 3 counts ~ outcome + treatment
  .. ..- attr(*, "variables")= language list(counts, outcome, treatment)
  .. ..- attr(*, "factors")= int [1:3, 1:2] 0 1 0 0 0 1
  .. .. ..- attr(*, "dimnames")=List of 2
  .. ..- attr(*, "term.labels")= chr [1:2] "outcome" "treatment"
  .. ..- attr(*, "order")= int [1:2] 1 1
  .. ..- attr(*, "intercept")= int 1
  .. ..- attr(*, "response")= int 1
  .. ..- attr(*, ".Environment")=<environment: R_GlobalEnv> 
  .. ..- attr(*, "predvars")= language list(counts, outcome, treatment)
  .. ..- attr(*, "dataClasses")= Named chr [1:3] "numeric" "factor" "factor"
  .. .. ..- attr(*, "names")= chr [1:3] "counts" "outcome" "treatment"
 $ family        :List of 12
  ..- attr(*, "class")= chr "family"
 $ deviance      : num 5.13
 $ aic           : num 56.8
 $ contrasts     :List of 2
 $ df.residual   : int 4
 $ null.deviance : num 10.6
 $ df.null       : int 8
 $ iter          : int 4
 $ deviance.resid: Named num [1:9] -0.671 0.963 -0.17 -0.22 -0.956 ...
  ..- attr(*, "names")= chr [1:9] "1" "2" "3" "4" ...
 $ coefficients  : num [1:5, 1:4] 3.04 -4.54e-01 -2.93e-01 1.34e-15 1.42e-15 ...
  ..- attr(*, "dimnames")=List of 2
 $ aliased       : Named logi [1:5] FALSE FALSE FALSE FALSE FALSE
  ..- attr(*, "names")= chr [1:5] "(Intercept)" "outcome2" "outcome3" "treatment2" ...
 $ dispersion    : num 1
 $ df            : int [1:3] 5 4 5
 $ cov.unscaled  : num [1:5, 1:5] 0.0292 -0.0159 -0.0159 -0.02 -0.02 ...
  ..- attr(*, "dimnames")=List of 2
 $ cov.scaled    : num [1:5, 1:5] 0.0292 -0.0159 -0.0159 -0.02 -0.02 ...
  ..- attr(*, "dimnames")=List of 2
 - attr(*, "class")= chr "summary.glm"

Por lo tanto, notamos el componente coefficients que podemos extraer usando coef(), pero otros componentes no tienen extractores, como null.deviance, que puede extraer como summ$null.deviance.