CODEWITHSUNDEEP
rglmp-value
ch-pub
ch-pub

Reputation: 1724

Extract pvalue from glm

I'm running many regressions and am only interested in the effect on the coefficient and p-value of one particular variable. So, in my script, I'd like to be able to just extract the p-value from the glm summary (getting the coefficient itself is easy). The only way I know of to view the p-value is using summary(myReg). Is there some other way?

e.g.:

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

I've tried treating fit$coefficients as a matrix, but am still unable to simply extract the p-value.

Is it possible to do this?

Thanks!

Upvotes: 53

Views: 89258

Answers (5)

user16247103
user16247103

Reputation: 31

The tidy function from the broom package (part of the Tidyverse, available on CRAN) provides a quick and easy way to convert your GLM summaries into a data frame, which might be useful in a number of situations other than the one you described above.

In this case, your desired output could be obtained with the code:

x1pValue <- broom::tidy(fit)$p.value[2]

Upvotes: 3

Ebby
Ebby

Reputation: 141

Well, this would be another way, however not the most efficient way to perform it :

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

This ensures that the extracted values from the glm is not in StatsBase. Therein, you can play around with pVals as per your requirement. Hope it helps, Ebby

Upvotes: 1

R. Prost
R. Prost

Reputation: 2088

Instead of the number you can put directly the name

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

the other ones available from the coefficient summary:

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

Upvotes: 12

Gavin Simpson
Gavin Simpson

Reputation: 174778

You want

coef(summary(fit))[,4]

which extracts the column vector of p values from the tabular output shown by summary(fit). The p-values aren't actually computed until you run summary() on the model fit.

By the way, use extractor functions rather than delve into objects if you can:

fit$coefficients[2]

should be

coef(fit)[2]

If there aren't extractor functions, str() is your friend. It allows you to look at the structure of any object, which allows you to see what the object contains and how to extract it:

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"

Hence we note the coefficients component which we can extract using coef(), but other components don't have extractors, like null.deviance, which you can extract as summ$null.deviance.

Upvotes: 76

James King
James King

Reputation: 6355

I've used this technique in the past to pull out predictor data from summary or from a fitted model object:

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

which lets me easily edit which variable I want to get data on.

Or as pointed out be @Ben, use match or %in%, somewhat cleaner than grepl:

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

Upvotes: 5

Related Questions