Data Science: предварительный анализ данных в R

Data Science: инструментарий и жизненный цикл проекта
Глушков Егор Александрович, гр. 20.М04-мм

Загрузка и описание данных

data(mtcars)
##    Переменная                                 Описание                  Шкала
## 1         mpg                              Миль/Галлон              Отношений
## 2         cyl                          Число цилиндров           Интервальная
## 3        disp              Рабочий объем (куб. дюймов)              Отношений
## 4          hp          Полная мощность (лошадиных сил)              Отношений
## 5        drat        Передаточное отношение задней оси              Отношений
## 6          wt                        Вес (1000 фунтов)              Отношений
## 7        qsec                           Время 1/4 мили              Отношений
## 8          vs       Двигатель (0=V-образный, 1=прямой) Номинальная (бинарная)
## 9          am Трансмиссия (0=автоматическая, 1=ручная) Номинальная (бинарная)
## 10       gear        Количество передач переднего хода           Интервальная
## 11       carb                  Количество карбюраторов           Интервальная
mtcars2 <- within(mtcars, {
   vs <- factor(vs, labels = c("V", "S"))
   am <- factor(am, labels = c("automatic", "manual"))
   cyl  <- ordered(cyl)
   gear <- ordered(gear)
   carb <- ordered(carb)
})
summary(mtcars2)
##       mpg        cyl         disp             hp             drat      
##  Min.   :10.40   4:11   Min.   : 71.1   Min.   : 52.0   Min.   :2.760  
##  1st Qu.:15.43   6: 7   1st Qu.:120.8   1st Qu.: 96.5   1st Qu.:3.080  
##  Median :19.20   8:14   Median :196.3   Median :123.0   Median :3.695  
##  Mean   :20.09          Mean   :230.7   Mean   :146.7   Mean   :3.597  
##  3rd Qu.:22.80          3rd Qu.:326.0   3rd Qu.:180.0   3rd Qu.:3.920  
##  Max.   :33.90          Max.   :472.0   Max.   :335.0   Max.   :4.930  
##        wt             qsec       vs             am     gear   carb  
##  Min.   :1.513   Min.   :14.50   V:18   automatic:19   3:15   1: 7  
##  1st Qu.:2.581   1st Qu.:16.89   S:14   manual   :13   4:12   2:10  
##  Median :3.325   Median :17.71                         5: 5   3: 3  
##  Mean   :3.217   Mean   :17.85                                4:10  
##  3rd Qu.:3.610   3rd Qu.:18.90                                6: 1  
##  Max.   :5.424   Max.   :22.90                                8: 1
head(mtcars2)
##                    mpg cyl disp  hp drat    wt  qsec vs        am gear carb
## Mazda RX4         21.0   6  160 110 3.90 2.620 16.46  V    manual    4    4
## Mazda RX4 Wag     21.0   6  160 110 3.90 2.875 17.02  V    manual    4    4
## Datsun 710        22.8   4  108  93 3.85 2.320 18.61  S    manual    4    1
## Hornet 4 Drive    21.4   6  258 110 3.08 3.215 19.44  S automatic    3    1
## Hornet Sportabout 18.7   8  360 175 3.15 3.440 17.02  V automatic    3    2
## Valiant           18.1   6  225 105 2.76 3.460 20.22  S automatic    3    1

Встроенный набор данных mtcars содержит \(32\) наблюдения по \(11\) признакам. Некоторые признаки перекодируем в соответствующие шкалы (интервальные, бинарные).

Проверка гипотез однородности

Проверим однородность подвыборок (признак Миль/Галлон mt) по фактору Тип двигателя (vs). Построим гистограмму частот с разбиением по указанному фактору.

library(ggplot2)

ggplot(data=mtcars2, aes(x=mpg, fill=vs)) +
  geom_histogram(binwidth = 1.25) +
  labs(title="Разбиение данных", 
       x="Миль/Галлон (mpg)", y="Количество") +
  scale_fill_discrete(name = "Двигатель (vs)", labels = c("V-образный", "Прямой (S)"))

Проиллюстрируем распределения в подвыборках с помощью “ящика с усами” (boxplot):

mtc <- data.frame(group=mtcars2$vs, X=mtcars2$mpg)
boxplot(X ~ group, mtc, xlab = "Тип двигателя", ylab = "Миль/Галлон")

Проверим нормальность данных в подвыборках:

p.Sh <- with(mtc, tapply(X, group, function(x)shapiro.test(x)$p.value)); p.Sh
##         V         S 
## 0.4491490 0.1665898

Согласие с нормальным распределением для обеих групп не отвергается в соответствии с критерием Шапиро-Уилка с \(p.value = 0.449\) и \(0.167\) при уровне значимости \(\alpha = 0.05\).

p.F <- var.test(X~group, mtc)$p.value; p.F
## [1] 0.1997209

В соответствии с критерием Фишера гипотеза о равенстве дисперсий двух групп не отвергается (\(p.value=0.1997\)), значит, можно использовать критерий Стьюдента (с параметром о равенстве дисперсий).

p.T <- t.test(X~group, mtc, var.equal=TRUE); p.T
## 
##  Two Sample t-test
## 
## data:  X by group
## t = -4.8644, df = 30, p-value = 3.416e-05
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##  -11.274221  -4.606732
## sample estimates:
## mean in group V mean in group S 
##        16.61667        24.55714

По критерию Стьюдента гипотеза о равенстве средних не принимается (\(p.value=3.415937 \times 10^{-5}\)) на уровне значимости \(\alpha = 0.05\), то есть различие средних в двух группах (различающихся по типу двигателя) значимо: автомобили с V-образным двигателем (\(16.62\) миль/галлон) действительно проезжают в среднем меньше, чем автомобили с прямым двигателем (\(24.56\)).

Множественная линейная регрессия

Построим модель множественной линейной регресии, зависящей от всех имеющихся переменных:

init.model <- lm(mpg~., data=mtcars)
summary(init.model)
## 
## Call:
## lm(formula = mpg ~ ., data = mtcars)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.4506 -1.6044 -0.1196  1.2193  4.6271 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)  
## (Intercept) 12.30337   18.71788   0.657   0.5181  
## cyl         -0.11144    1.04502  -0.107   0.9161  
## disp         0.01334    0.01786   0.747   0.4635  
## hp          -0.02148    0.02177  -0.987   0.3350  
## drat         0.78711    1.63537   0.481   0.6353  
## wt          -3.71530    1.89441  -1.961   0.0633 .
## qsec         0.82104    0.73084   1.123   0.2739  
## vs           0.31776    2.10451   0.151   0.8814  
## am           2.52023    2.05665   1.225   0.2340  
## gear         0.65541    1.49326   0.439   0.6652  
## carb        -0.19942    0.82875  -0.241   0.8122  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.65 on 21 degrees of freedom
## Multiple R-squared:  0.869,  Adjusted R-squared:  0.8066 
## F-statistic: 13.93 on 10 and 21 DF,  p-value: 3.793e-07

Согласно \(p.value=3.793 \times 10^{-7}\) с уровнем значимости \(\alpha = 0.001\) не принимается гипотеза о том, что данная модель не зависит от переменных (предикторов). При этом лишь переменная wt отлична от нуля с уровнем значимости \(0.1\). Коэффициент детерминации равен \(0.869\), скорректированный – \(0.8066\).

Произведем пошаговую подгонку регрессии, когда на каждом шаге добавляется наиболее значимая переменная (в соответствии с информационным критерием Акаике AIC).

library(MASS)

empty_m <- lm(mpg ~ 1, mtcars)
SA.f <- stepAIC(empty_m, direction = 'forward', scope=list(lower=empty_m, upper=~cyl + disp + hp + drat + wt + qsec + vs + am + gear + carb))
## Start:  AIC=115.94
## mpg ~ 1
## 
##        Df Sum of Sq     RSS     AIC
## + wt    1    847.73  278.32  73.217
## + cyl   1    817.71  308.33  76.494
## + disp  1    808.89  317.16  77.397
## + hp    1    678.37  447.67  88.427
## + drat  1    522.48  603.57  97.988
## + vs    1    496.53  629.52  99.335
## + am    1    405.15  720.90 103.672
## + carb  1    341.78  784.27 106.369
## + gear  1    259.75  866.30 109.552
## + qsec  1    197.39  928.66 111.776
## <none>              1126.05 115.943
## 
## Step:  AIC=73.22
## mpg ~ wt
## 
##        Df Sum of Sq    RSS    AIC
## + cyl   1    87.150 191.17 63.198
## + hp    1    83.274 195.05 63.840
## + qsec  1    82.858 195.46 63.908
## + vs    1    54.228 224.09 68.283
## + carb  1    44.602 233.72 69.628
## + disp  1    31.639 246.68 71.356
## <none>              278.32 73.217
## + drat  1     9.081 269.24 74.156
## + gear  1     1.137 277.19 75.086
## + am    1     0.002 278.32 75.217
## 
## Step:  AIC=63.2
## mpg ~ wt + cyl
## 
##        Df Sum of Sq    RSS    AIC
## + hp    1   14.5514 176.62 62.665
## + carb  1   13.7724 177.40 62.805
## <none>              191.17 63.198
## + qsec  1   10.5674 180.60 63.378
## + gear  1    3.0281 188.14 64.687
## + disp  1    2.6796 188.49 64.746
## + vs    1    0.7059 190.47 65.080
## + am    1    0.1249 191.05 65.177
## + drat  1    0.0010 191.17 65.198
## 
## Step:  AIC=62.66
## mpg ~ wt + cyl + hp
## 
##        Df Sum of Sq    RSS    AIC
## <none>              176.62 62.665
## + am    1    6.6228 170.00 63.442
## + disp  1    6.1762 170.44 63.526
## + carb  1    2.5187 174.10 64.205
## + drat  1    2.2453 174.38 64.255
## + qsec  1    1.4010 175.22 64.410
## + gear  1    0.8558 175.76 64.509
## + vs    1    0.0599 176.56 64.654
summary(SA.f)
## 
## Call:
## lm(formula = mpg ~ wt + cyl + hp, data = mtcars)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.9290 -1.5598 -0.5311  1.1850  5.8986 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 38.75179    1.78686  21.687  < 2e-16 ***
## wt          -3.16697    0.74058  -4.276 0.000199 ***
## cyl         -0.94162    0.55092  -1.709 0.098480 .  
## hp          -0.01804    0.01188  -1.519 0.140015    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.512 on 28 degrees of freedom
## Multiple R-squared:  0.8431, Adjusted R-squared:  0.8263 
## F-statistic: 50.17 on 3 and 28 DF,  p-value: 2.184e-11

Итоговая модель построена на трёх переменных: wt, cyl, hp, первые две из которых значимо отличны от нуля (с уровнями значимости \(0.001\), \(0.1\) соответственно). Сама модель также значимо зависит от переменных с \(p.value=2.18 \times 10^{-11}\).

Коэффициент детерминации равен \(0.8431\), скорректированный – \(0.8263\).

Модель с “обратным ходом”, исключающая переменные:

SA.b <- stepAIC(init.model, direction = 'backward')
## Start:  AIC=70.9
## mpg ~ cyl + disp + hp + drat + wt + qsec + vs + am + gear + carb
## 
##        Df Sum of Sq    RSS    AIC
## - cyl   1    0.0799 147.57 68.915
## - vs    1    0.1601 147.66 68.932
## - carb  1    0.4067 147.90 68.986
## - gear  1    1.3531 148.85 69.190
## - drat  1    1.6270 149.12 69.249
## - disp  1    3.9167 151.41 69.736
## - hp    1    6.8399 154.33 70.348
## - qsec  1    8.8641 156.36 70.765
## <none>              147.49 70.898
## - am    1   10.5467 158.04 71.108
## - wt    1   27.0144 174.51 74.280
## 
## Step:  AIC=68.92
## mpg ~ disp + hp + drat + wt + qsec + vs + am + gear + carb
## 
##        Df Sum of Sq    RSS    AIC
## - vs    1    0.2685 147.84 66.973
## - carb  1    0.5201 148.09 67.028
## - gear  1    1.8211 149.40 67.308
## - drat  1    1.9826 149.56 67.342
## - disp  1    3.9009 151.47 67.750
## - hp    1    7.3632 154.94 68.473
## <none>              147.57 68.915
## - qsec  1   10.0933 157.67 69.032
## - am    1   11.8359 159.41 69.384
## - wt    1   27.0280 174.60 72.297
## 
## Step:  AIC=66.97
## mpg ~ disp + hp + drat + wt + qsec + am + gear + carb
## 
##        Df Sum of Sq    RSS    AIC
## - carb  1    0.6855 148.53 65.121
## - gear  1    2.1437 149.99 65.434
## - drat  1    2.2139 150.06 65.449
## - disp  1    3.6467 151.49 65.753
## - hp    1    7.1060 154.95 66.475
## <none>              147.84 66.973
## - am    1   11.5694 159.41 67.384
## - qsec  1   15.6830 163.53 68.200
## - wt    1   27.3799 175.22 70.410
## 
## Step:  AIC=65.12
## mpg ~ disp + hp + drat + wt + qsec + am + gear
## 
##        Df Sum of Sq    RSS    AIC
## - gear  1     1.565 150.09 63.457
## - drat  1     1.932 150.46 63.535
## <none>              148.53 65.121
## - disp  1    10.110 158.64 65.229
## - am    1    12.323 160.85 65.672
## - hp    1    14.826 163.35 66.166
## - qsec  1    26.408 174.94 68.358
## - wt    1    69.127 217.66 75.350
## 
## Step:  AIC=63.46
## mpg ~ disp + hp + drat + wt + qsec + am
## 
##        Df Sum of Sq    RSS    AIC
## - drat  1     3.345 153.44 62.162
## - disp  1     8.545 158.64 63.229
## <none>              150.09 63.457
## - hp    1    13.285 163.38 64.171
## - am    1    20.036 170.13 65.466
## - qsec  1    25.574 175.67 66.491
## - wt    1    67.572 217.66 73.351
## 
## Step:  AIC=62.16
## mpg ~ disp + hp + wt + qsec + am
## 
##        Df Sum of Sq    RSS    AIC
## - disp  1     6.629 160.07 61.515
## <none>              153.44 62.162
## - hp    1    12.572 166.01 62.682
## - qsec  1    26.470 179.91 65.255
## - am    1    32.198 185.63 66.258
## - wt    1    69.043 222.48 72.051
## 
## Step:  AIC=61.52
## mpg ~ hp + wt + qsec + am
## 
##        Df Sum of Sq    RSS    AIC
## - hp    1     9.219 169.29 61.307
## <none>              160.07 61.515
## - qsec  1    20.225 180.29 63.323
## - am    1    25.993 186.06 64.331
## - wt    1    78.494 238.56 72.284
## 
## Step:  AIC=61.31
## mpg ~ wt + qsec + am
## 
##        Df Sum of Sq    RSS    AIC
## <none>              169.29 61.307
## - am    1    26.178 195.46 63.908
## - qsec  1   109.034 278.32 75.217
## - wt    1   183.347 352.63 82.790
summary(SA.b)
## 
## Call:
## lm(formula = mpg ~ wt + qsec + am, data = mtcars)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.4811 -1.5555 -0.7257  1.4110  4.6610 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   9.6178     6.9596   1.382 0.177915    
## wt           -3.9165     0.7112  -5.507 6.95e-06 ***
## qsec          1.2259     0.2887   4.247 0.000216 ***
## am            2.9358     1.4109   2.081 0.046716 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.459 on 28 degrees of freedom
## Multiple R-squared:  0.8497, Adjusted R-squared:  0.8336 
## F-statistic: 52.75 on 3 and 28 DF,  p-value: 1.21e-11

Модель с “обратным ходом” значимо (\(p.value=1.21 \times 10^{-11}\)) зависит от трех переменных (wt, qsec, am), каждая из которых является ненулевой при уровне значимости \(\alpha = 0.05\). Коэффициент детерминации равен \(0.8497\), скорректированный – \(0.8336\).

Таким образом, последняя модель множественной линейной регрессии, исключающая переменные, показывает наилучший результат (в соответствии со скорректированным коэффициентом детерминации).