Statistisko modeļu rezultātu pārvēršana par datu tabulām ar paketi broom
Problēma
Veicot statistiskās analīzes programmā R, piemēram, lineāro regresiju vai T testu, mēs iegūstam ļoti saprotamu rezultātu, kas parāda mūs interesējošās lietas.
Šeit ir piemērs regresijas analīzei izmantojot iris datu objektu starp mainīgajiem Petal.Width un Sepal.Width. Analīzes rezultātos mēs redzam gan aprēķinātos koeficientus, gan arī to būtiskumu.
data(iris)
mod<-lm(Petal.Width ~ Sepal.Width, data = iris)
summary(mod)##
## Call:
## lm(formula = Petal.Width ~ Sepal.Width, data = iris)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.38424 -0.60889 -0.03208 0.52691 1.64812
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.1569 0.4131 7.642 2.47e-12 ***
## Sepal.Width -0.6403 0.1338 -4.786 4.07e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.7117 on 148 degrees of freedom
## Multiple R-squared: 0.134, Adjusted R-squared: 0.1282
## F-statistic: 22.91 on 1 and 148 DF, p-value: 4.073e-06Problēmas rodas tajā brīdī, kad mums šos rezultātus vajag dabūt kā datu tabulu, lai eksportētu vai arī izmantotu kādās citās analīzes/attēlu veidošanā. Vairumā gadījumu statistisko testu rezultāti ir saglabāti kā saraksti (list) ar ļoti daudz apakšelementiem, kas apgrūtina rezultātu tālāku izmantošanu.
Risinājums
broom pakete
Ar šo problēmu veiksmīgi var tikt galā izmantojot paketi broom, kuras trīs galvenās funkcijas ir tidy(), augment() un glance().
tidy()
Funkcija tidy() ir paredzēta, lai kā datu tabulu izveidotu koeficientus un to būtiskumu, ko iegūstam, piemēram, lineārā regresijas rezultātā. Tagad katrs koeficients ir savā rindā (kolonna term), kā arī ir kolonnas ar pašiem koeficientiem, standartkļūdām un p vērtībām. Turklāt kolonnu nosaukumi nesatur dīvainus apzīmējumus vai atstarpes.
library(broom)
tidy(lm(Petal.Width ~ Sepal.Width, data = iris))## # A tibble: 2 x 5
## term estimate std.error statistic p.value
## <chr> <dbl> <dbl> <dbl> <dbl>
## 1 (Intercept) 3.16 0.413 7.64 2.47e-12
## 2 Sepal.Width -0.640 0.134 -4.79 4.07e- 6Līdzīgi mēs varam darīt arī ar T testa rezultātiem.
tidy(t.test(Sepal.Width~Species,data=subset(iris,Species!="setosa")))## # A tibble: 1 x 10
## estimate estimate1 estimate2 statistic p.value parameter conf.low conf.high
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 -0.204 2.77 2.97 -3.21 0.00182 97.9 -0.330 -0.0777
## # … with 2 more variables: method <chr>, alternative <chr>glance()
Ja ir nepieciešams kā datu tabulu saglabāt dažādus modeļa kopējos novērtējuma radītājus, piemēram, R2, p vērtību, AIC, tad ir jāizmanto funkcija glance().
glance(lm(Petal.Width ~ Sepal.Width, data = iris))## # A tibble: 1 x 11
## r.squared adj.r.squared sigma statistic p.value df logLik AIC BIC
## <dbl> <dbl> <dbl> <dbl> <dbl> <int> <dbl> <dbl> <dbl>
## 1 0.134 0.128 0.712 22.9 4.07e-6 2 -161. 328. 337.
## # … with 2 more variables: deviance <dbl>, df.residual <int>augment()
Ar funkcijas augment() palīdzību ir iespējams no modeļiem izvilkt katram novērojumam aprēķināmos lielumus un tos pievienot modeļa veidošanai izmantotajiem datiem. Tas ietver prognozētās vērtības, atlikuma vērtības, kā arī diagnosticējošos rādītājus.
head(augment(lm(Petal.Width ~ Sepal.Width, data = iris)))## # A tibble: 6 x 9
## Petal.Width Sepal.Width .fitted .se.fit .resid .hat .sigma .cooksd
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 0.2 3.5 0.916 0.0830 -0.716 0.0136 0.712 0.00707
## 2 0.2 3 1.24 0.0586 -1.04 0.00678 0.709 0.00729
## 3 0.2 3.2 1.11 0.0612 -0.908 0.00739 0.710 0.00610
## 4 0.2 3.1 1.17 0.0584 -0.972 0.00673 0.710 0.00636
## 5 0.2 3.6 0.852 0.0930 -0.652 0.0171 0.712 0.00741
## 6 0.4 3.9 0.660 0.127 -0.260 0.0318 0.714 0.00226
## # … with 1 more variable: .std.resid <dbl>Saites
Lai iegūtu papildus informāciju par šo paketi un tās darbību, var apmeklēt:
Didzis Elferts
Datu zinātnieks, vadošais pētnieks
Manas zinātniskās intereses saistītas ar statistiku un programmu R.