what is the code to extract SE in regression?

chiaoyi · 2020-02-13T02:44:17.169Z

I know coef is for coefficients. Now I want to extract standard error.

nirgrahamuk · 2020-02-13T02:45:45.266Z

If you have an object you can see what's in it with str() function. It can be useful to see how things are named. Does that help in this case?

chiaoyi · 2020-02-13T02:48:28.053Z

I'm not sure which one is.

nirgrahamuk · 2020-02-13T02:50:50.674Z

What object did you put your regression into?
or
What code do you use to get your coef?

chiaoyi · 2020-02-13T02:52:45.207Z

coef(reg.ldv)["treatment"]
Now I want the SE of treatment variable.

nirgrahamuk · 2020-02-13T02:55:31.737Z

Please. What is result of

str(reg.ldv)

williaml · 2020-02-13T02:56:48.708Z

You can also use the broom package to get the results of your regression into a data frame.

for example:

glmfit <- glm(am ~ wt, mtcars, family="binomial")
tidy(glmfit)

## # A tibble: 2 x 5
##   term        estimate std.error statistic p.value
##   <chr>          <dbl>     <dbl>     <dbl>   <dbl>
## 1 (Intercept)    12.0       4.51      2.67 0.00759
## 2 wt             -4.02      1.44     -2.80 0.00509

https://cran.r-project.org/web/packages/broom/vignettes/broom.html

chiaoyi · 2020-02-13T02:57:17.236Z

> str(reg.ldv)
List of 13
 $ coefficients : Named num [1:3] 7.75 -0.362 0.891
  ..- attr(*, "names")= chr [1:3] "(Intercept)" "treatment" "daily_kwh_lag"
 $ residuals    : Named num [1:452155] 3.66 -1.71 -7.17 2.79 -5.03 ...
  ..- attr(*, "names")= chr [1:452155] "1" "2" "3" "4" ...
 $ effects      : Named num [1:452155] -37588.71 137.32 14590.08 2.79 -5.02 ...
  ..- attr(*, "names")= chr [1:452155] "(Intercept)" "treatment" "daily_kwh_lag" "" ...
 $ rank         : int 3
 $ fitted.values: Named num [1:452155] 51.9 52.8 51.3 29.1 38.4 ...
  ..- attr(*, "names")= chr [1:452155] "1" "2" "3" "4" ...
 $ assign       : int [1:3] 0 1 2
 $ qr           :List of 5
  ..$ qr   : num [1:452155, 1:3] -672.42472 0.00149 0.00149 0.00149 0.00149 ...
  .. ..- attr(*, "dimnames")=List of 2
  .. .. ..$ : chr [1:452155] "1" "2" "3" "4" ...
  .. .. ..$ : chr [1:3] "(Intercept)" "treatment" "daily_kwh_lag"
  .. ..- attr(*, "assign")= int [1:3] 0 1 2
  ..$ qraux: num [1:3] 1 1 1
  ..$ pivot: int [1:3] 1 2 3
  ..$ tol  : num 0.0000001
  ..$ rank : int 3
  ..- attr(*, "class")= chr "qr"
 $ df.residual  : int 452152
 $ na.action    : 'omit' Named int [1:145] 3859 3860 7899 18506 21906 22785 27730 30457 40087 40873 ...
  ..- attr(*, "names")= chr [1:145] "3859" "3860" "7899" "18506" ...
 $ xlevels      : Named list()
 $ call         : language lm(formula = daily_kwh ~ treatment + daily_kwh_lag, data = df.post)
 $ terms        :Classes 'terms', 'formula'  language daily_kwh ~ treatment + daily_kwh_lag
  .. ..- attr(*, "variables")= language list(daily_kwh, treatment, daily_kwh_lag)
  .. ..- attr(*, "factors")= int [1:3, 1:2] 0 1 0 0 0 1
  .. .. ..- attr(*, "dimnames")=List of 2
  .. .. .. ..$ : chr [1:3] "daily_kwh" "treatment" "daily_kwh_lag"
  .. .. .. ..$ : chr [1:2] "treatment" "daily_kwh_lag"
  .. ..- attr(*, "term.labels")= chr [1:2] "treatment" "daily_kwh_lag"
  .. ..- attr(*, "order")= int [1:2] 1 1
  .. ..- attr(*, "intercept")= int 1
  .. ..- attr(*, "response")= int 1
  .. ..- attr(*, ".Environment")=<environment: R_GlobalEnv> 
  .. ..- attr(*, "predvars")= language list(daily_kwh, treatment, daily_kwh_lag)
  .. ..- attr(*, "dataClasses")= Named chr [1:3] "numeric" "numeric" "numeric"
  .. .. ..- attr(*, "names")= chr [1:3] "daily_kwh" "treatment" "daily_kwh_lag"
 $ model        :'data.frame':	452155 obs. of  3 variables:
  ..$ daily_kwh    : num [1:452155] 55.6 51.1 44.1 31.9 33.3 ...
  ..$ treatment    : int [1:452155] 1 1 1 1 1 1 1 1 1 1 ...
  ..$ daily_kwh_lag: num [1:452155] 50 51 49.2 24.3 34.8 ...
  ..- attr(*, "terms")=Classes 'terms', 'formula'  language daily_kwh ~ treatment + daily_kwh_lag
  .. .. ..- attr(*, "variables")= language list(daily_kwh, treatment, daily_kwh_lag)
  .. .. ..- attr(*, "factors")= int [1:3, 1:2] 0 1 0 0 0 1
  .. .. .. ..- attr(*, "dimnames")=List of 2
  .. .. .. .. ..$ : chr [1:3] "daily_kwh" "treatment" "daily_kwh_lag"
  .. .. .. .. ..$ : chr [1:2] "treatment" "daily_kwh_lag"
  .. .. ..- attr(*, "term.labels")= chr [1:2] "treatment" "daily_kwh_lag"
  .. .. ..- attr(*, "order")= int [1:2] 1 1
  .. .. ..- attr(*, "intercept")= int 1
  .. .. ..- attr(*, "response")= int 1
  .. .. ..- attr(*, ".Environment")=<environment: R_GlobalEnv> 
  .. .. ..- attr(*, "predvars")= language list(daily_kwh, treatment, daily_kwh_lag)
  .. .. ..- attr(*, "dataClasses")= Named chr [1:3] "numeric" "numeric" "numeric"
  .. .. .. ..- attr(*, "names")= chr [1:3] "daily_kwh" "treatment" "daily_kwh_lag"
  ..- attr(*, "na.action")= 'omit' Named int [1:145] 3859 3860 7899 18506 21906 22785 27730 30457 40087 40873 ...
  .. ..- attr(*, "names")= chr [1:145] "3859" "3860" "7899" "18506" ...
 - attr(*, "class")= chr "lm"

nirgrahamuk · 2020-02-13T03:01:23.901Z

Williaml gave good solution if the model type is suitable to what you are doing.
Can try summary() on your model. To see what you get maybe.
It might be that something like this possible:

#Residual Standard error (Like Standard Deviation)
k=length(model$coefficients)-1 #Subtract one to ignore intercept
SSE=sum(model$residuals**2)
n=length(model$residuals)
sqrt(SSE/(n-(1+k))) #Residual Standard Error

chiaoyi · 2020-02-13T04:40:55.132Z

no single code to extract standard error?

williaml · 2020-02-13T05:12:39.915Z

What's wrong with extracting it from the dataframe?

For example:

library(broom)
library(tidyverse)
fit <- lm(am ~ wt, mtcars)
a <- tidy(fit)
a %>% 
  select(std.error) %>% 
  pull() # just the values in a vector

[1] 0.22557907 0.06717067

system · 2020-02-20T05:12:40.527Z

This topic was automatically closed 7 days after the last reply. New replies are no longer allowed.