28 февраля 2012

Расчет параметров описательной статистики в R




Использование специальных функций

Благодаря наличию специально созданных для этого функций, расчет параметров описательной статистики в R не составляет никакого труда. Ниже я продемонстрирую использование этих функций на примере ранее рассмотренных данных по характеристикам 32 моделей автомобилей (таблица mtcars, входящая в стандартный набор данных R):

data(mtcars)
mtcars
                     mpg cyl  disp  hp drat    wt  qsec vs am gear carb
Mazda RX4           21.0   6 160.0 110 3.90 2.620 16.46  0  1    4    4
Mazda RX4 Wag       21.0   6 160.0 110 3.90 2.875 17.02  0  1    4    4
Datsun 710          22.8   4 108.0  93 3.85 2.320 18.61  1  1    4    1
Hornet 4 Drive      21.4   6 258.0 110 3.08 3.215 19.44  1  0    3    1
Hornet Sportabout   18.7   8 360.0 175 3.15 3.440 17.02  0  0    3    2
Valiant             18.1   6 225.0 105 2.76 3.460 20.22  1  0    3    1
Duster 360          14.3   8 360.0 245 3.21 3.570 15.84  0  0    3    4
Merc 240D           24.4   4 146.7  62 3.69 3.190 20.00  1  0    4    2
Merc 230            22.8   4 140.8  95 3.92 3.150 22.90  1  0    4    2
Merc 280            19.2   6 167.6 123 3.92 3.440 18.30  1  0    4    4
...
...


Каждая модель описана по 11 признакам, два из которых (vs и am) являются номинальными переменными (факторами) c уровнями, закодированными в виде 0 и 1 (подробнее см. ?mtcars).

Для расчета  арифметической средней, медианы, дисперсии, стандартного отклонения, а также минимального и максимального значений в R служат функции mean(), median(), var(), sd(), min() и max() соответственно. Используем эти функции в отношении, например, параметра mpg (пробег автомобиля (в милях) в расчете на один галлон топлива):

# Арифметическая средняя:
mean(mtcars$mpg)
 
[1] 20.1
 
# Медиана
median(mtcars$mpg)
[1] 19.2
 
# Дисперсия:
var(mtcars$mpg)
[1] 36.3
 
# Стандартное отклонение:
sd(mtcars$mpg)
[1] 6.0
# Минимальное значение:
min(mtcars$mpg)
[1] 10.4
 
# Максимальное значение:
max(mtcars$mpg)
[1] 33.9

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



На языке R мы можем записать это следующим образом:

SEmpg = sd(mtcars$mpg)/sqrt(length(mtcars$mpg))
# функция length() возвращает число элементов в векторе mpg

Получаем:

SEmpg
[1] 1.065

Квантили рассчитываются в R при помощи функции quantile():

quantile(mtcars$mpg)
    0%    25%    50%    75%   100% 
10.400 15.425 19.200 22.800 33.900

При настройках, заданных по умолчанию, выполнение указанной команды приведет к расчету минимального (10.4) и максимального (33.9) значений, а также трех квартилей, т.е. значений, которые делят совокупность на четыре равные части - 15.4, 19.2 и 22.8. Разница между первым и третим квартилями носит название интерквартильный размах (ИКР; англ. interquartile range). ИКР является робастным аналогом дисперсии и может быть рассчитан в R при помощи функции IQR():

IQR(mtcars$mpg)
[1] 7.375

Функция quantile() позволяет рассчитать и другие квантили. Например, децили (т.е. значения, делящие совокупность на десять частей) можно получить следующим образом:

quantile(mtcars$mpg, p = seq(0, 1, 0.1))
   0%   10%   20%   30%   40%   50%   60%   70%   80%   90%  100% 
10.40 14.34 15.20 15.98 17.92 19.20 21.00 21.47 24.08 30.09 33.90

В приведенной команде важен аргумент p (от probability - вероятность), при помощи которого был задан вектор чисел от 0 до 1 с шагом 0.1.

Обратите внимание на то, что существует несколько способов оценки квантилей по выборочным данным. Подробнее об этом можно узнать в справочном файле по функции quantile() (доступен по команде ?quantile).

Отсутствующие значения в данных могут несколько усложнить вычисления. В качестве демонстрации заменим 3-е значение переменной mpg на NA (от not available - не доступно) - обозначение, используемое в R для отсутствующих наблюдений, - а затем попытаемся вычислить среднее значение:

mtcars$mpg[3] <- NA
 
# Просмотрим результат:
head(mtcars$mpg)
[1] 21.0 21.0   NA 21.4 18.7 18.1
 
# Попробуем рассчитать среднее значение для mpg:
mean(mtcars$mpg)
[1] NA

Ничего не вышло - вместо среднего значения программа выдала NA, что вполне логично. R не будет пропускать отсутствующие значения автоматически, если мы не включим соответствующую опцию - na.rm (от not available и remove - удалить):

mean(mtcars$mpg, na.rm = TRUE)
[1] 20.0

Рассмотренный прием срабатывает в большинстве случаев. Одним из немногих исключений является функция length(), используемая для определения размера вектора. Аргумент na.rm у этой функции отсутствует, так что подсчитаны будут и имеющиеся, и отсутствующие значения:

length(mtcars$mpg)
[1] 32

Если все же стоит задача подсчитать количество неотсутствующих значений, то можно воспользоваться следующим приемом:

sum(!is.na(mtcars$mpg))
[1] 31

Ключом здесь является использование команды is.na(mtcars$mpg), которая проверяет каждое значение mpg и возвращает FALSE, если это значение не равно NA, и TRUE иначе. В сочетании с логическим оператором ! ("не"), команда sum далее подсчитывает только те значения mpg, которые не равны NA (логические TRUE здесь конвертируются в 1, которые можно суммировать).

Существуют еще две функции, которые могут оказаться полезными при анализе свойств совокупностей - which.min() и which.max(). Как следует из названий, эти функции позволяют выяснить порядковые номера элементов, обладающих минимальным и максимальным значениями соответственно. Если минимальное/максимальное значение принимают несколько элементов в векторе, то будет возвращен порядковый номер первого элемента с этим значением. В случае с mpg имеем:

which.min(mtcars$mpg)
[1] 15
 
which.max(mtcars$mpg)
[1] 20

Видим, что минимальный и максимальный пробег в расчете на галлон топлива имеют модели под номерами 15 и 20 соответственно. Выяснить названия этих моделей мы можем, совместив команды which.min() и which.max() с командой rownames() (от row - строка, и names - имена):

rownames(mtcars)[which.min(mtcars$mpg)]
[1] "Cadillac Fleetwood"
 
rownames(mtcars)[which.max(mtcars$mpg)]
[1] "Toyota Corolla"

Подробнее о приемах индексирования векторов в R см. здесь.


Использование функции summary()

В системе R имеется возможность и более быстрого расчета основных параметров описательной статистики. Для этого, в частности, служит функция общего назначения summary():

summary(mtcars$mpg)
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    NA's 
  10.40   15.35   19.20   20.00   22.15   33.90    1.00

Всего одной строки кода оказалось достаточно для получения минимального (Min) и максимального (Max) значений переменной mpg, медианы (Median), арифметической средней (Mean), первого (1st Qu.) и третьего (3rd Qu.) квартилей, а также для выяснения количества отсутствующих значений (NA's). Более того, подобную сводку мы можем получить сразу для всей таблицы данных:

summary(mtcars)
      mpg             cyl             disp             hp       
 Min.   :10.40   Min.   :4.000   Min.   : 71.1   Min.   : 52.0  
 1st Qu.:15.35   1st Qu.:4.000   1st Qu.:120.8   1st Qu.: 96.5  
 Median :19.20   Median :6.000   Median :196.3   Median :123.0  
 Mean   :20.00   Mean   :6.188   Mean   :230.7   Mean   :146.7  
 3rd Qu.:22.15   3rd Qu.:8.000   3rd Qu.:326.0   3rd Qu.:180.0  
 Max.   :33.90   Max.   :8.000   Max.   :472.0   Max.   :335.0  
 NA's   : 1.00                                                  
      drat             wt             qsec             vs        
 Min.   :2.760   Min.   :1.513   Min.   :14.50   Min.   :0.0000  
 1st Qu.:3.080   1st Qu.:2.581   1st Qu.:16.89   1st Qu.:0.0000  
 Median :3.695   Median :3.325   Median :17.71   Median :0.0000  
 Mean   :3.597   Mean   :3.217   Mean   :17.85   Mean   :0.4375  
 3rd Qu.:3.920   3rd Qu.:3.610   3rd Qu.:18.90   3rd Qu.:1.0000  
 Max.   :4.930   Max.   :5.424   Max.   :22.90   Max.   :1.0000  
 
       am              gear            carb      
 Min.   :0.0000   Min.   :3.000   Min.   :1.000  
 1st Qu.:0.0000   1st Qu.:3.000   1st Qu.:2.000  
 Median :0.0000   Median :4.000   Median :2.000  
 Mean   :0.4062   Mean   :3.688   Mean   :2.812  
 3rd Qu.:1.0000   3rd Qu.:4.000   3rd Qu.:4.000  
 Max.   :1.0000   Max.   :5.000   Max.   :8.000 

Результат выглядел бы замечательно, если бы не одно "но". Переменные vs и am являются факторами, но уровни их закодированы при помощи чисел 0 и 1. К сожалению, система R не распознала эти две переменные как факторы и рассчитала соответствующие параметры описательной статистики, как для обычных числовых переменных. Однако мы можем изменить такое поведение R, самостоятельно преобразовав vs и am в факторы при помощи функции as.factor():

mtcars$vs <- as.factor(mtcars$vs)
 
mtcars$am <- as.factor(mtcars$am)
 
# Проверим, удалась ли конвертация:
is.factor(mtcars$vs)
[1] TRUE
 
is.factor(mtcars$am)
[1] TRUE

Теперь результат действия функции summary() в отношении таблицы mtcars будет выглядеть так:

summary(mtcars)
      mpg             cyl             disp             hp       
 Min.   :10.40   Min.   :4.000   Min.   : 71.1   Min.   : 52.0  
 1st Qu.:15.35   1st Qu.:4.000   1st Qu.:120.8   1st Qu.: 96.5  
 Median :19.20   Median :6.000   Median :196.3   Median :123.0  
 Mean   :20.00   Mean   :6.188   Mean   :230.7   Mean   :146.7  
 3rd Qu.:22.15   3rd Qu.:8.000   3rd Qu.:326.0   3rd Qu.:180.0  
 Max.   :33.90   Max.   :8.000   Max.   :472.0   Max.   :335.0  
 NA's   : 1.00
      drat             wt             qsec       vs     am    
 Min.   :2.760   Min.   :1.513   Min.   :14.50   0:18   0:19  
 1st Qu.:3.080   1st Qu.:2.581   1st Qu.:16.89   1:14   1:13  
 Median :3.695   Median :3.325   Median :17.71                
 Mean   :3.597   Mean   :3.217   Mean   :17.85                
 3rd Qu.:3.920   3rd Qu.:3.610   3rd Qu.:18.90                
 Max.   :4.930   Max.   :5.424   Max.   :22.90                
 
      gear            carb      
 Min.   :3.000   Min.   :1.000  
 1st Qu.:3.000   1st Qu.:2.000  
 Median :4.000   Median :2.000  
 Mean   :3.688   Mean   :2.812

Обратите внимание на сводки по vs и am: поскольку эти переменные теперь распознаны программой как факторы, единственный способ описать их - это подсчитать количество наблюдений для каждого уровня.


Использование функции tapply()

Функция tapply() принадлежит к важному "apply-семейству" R-функций. Эти функции позволяют выполнять математические вычисления над определенными элементами таблиц данных, матриц, или массивов (например, быстро вычислять среднее значение для каждого столбца или строки таблицы, и т.п.).

Предположим, мы хотим выяснить средний объем двигателя (переменная disp, в кубических дюймах) у моделей с автоматической и ручной коробкой передач (переменная am; 1 - ручная коробка, 0 - автоматическая коробка). Функция tapply() позволяет сделать это следующим образом:

tapply(X = mtcars$disp, INDEX = mtcars$am, FUN = mean)
     0          1 
290.38     143.53

Как видно из приведенной команды, основными аргументами функции tapply() являются:
  • Х - числовой вектор
  • INDEX - список факторов, для уровней которых рассчитываются значения функции FUN
  • FUN - любая, в том числе пользовательская, функция
Поскольку аргумент INDEX способен принимать список из нескольких факторов, мы можем усложнить приведенную выше команду:

tapply(X = mtcars$disp, INDEX = list(mtcars$am, mtcars$vs), FUN = mean)
       0      1
0 357.62 175.11
1 206.22  89.80

Аргумент FUN, как уже было отмечено, может принимать любые, в том числе и пользовательские функции. Рассмотрим, например, расчет стандартных ошибок для средних значений объема двигателя у автомобилей с автоматической и ручной коробкой передач. Для начала создадим функцию SE для расчета стандартных ошибок (см. пример выше):

SE <- function(x) {sd(x)/sqrt(length(x))}

Теперь совместим эту новую функцию с tapply():

tapply(X = mtcars$disp, INDEX = mtcars$am, FUN = SE)
   0    1 
25.3 24.2

Таким образом, объем двигателя у моделей с автоматической коробкой передач составляет в среднем 290.4 ± 25.3, а у автомобилей с механической коробкой - 143.5 ± 24.2 кубических дюймов.


Использование возможностей дополнительных пакетов

Рассмотренные выше функции позволяют получить достаточно полное представление об анализируемых выборках и таблицах данных. Однако специальные функции для расчета некоторых параметров описательной статистики не входят в базовую версию R. С одним из таких параметров мы уже столкнулись - стандартная ошибка арифметической средней. Другие примеры включают коэффициенты эксцесса (англ. kurtosis) и асимметрии (skewness) - параметры, характеризующие форму распределения. Конечно, мы можем рассчитать эти величины по соответствующим формулам или даже написать собственные функции для этих целей. Однако это уже было сделано до нас - достаточно воспользоваться имеющимися дополнительными пакетами для R, например, пакетом moments. Если этот пакет не установлен на Вашем компьютере, выполните следующую команду (естественно, Ваш компьютер должен быть при этом подключен к Internet):

install.packages("moments")

Рассчитать коэффициенты эксцесса и асимметрии теперь очень просто:

library(moments)  #загрузка пакета moments
 
kurtosis(mtcars$mpg, na.rm = TRUE)
[1] 2.79
 
skewness(mtcars$mpg, na.rm = TRUE)
[1] 0.68

Многие R-пакеты имеют собственные функции, аналогичные стандартной summary(), для вывода компактных описательных сводок по таблицам данных. Ниже приведены несколько примеров таких пакетов и функций (предполагается, что соответствующий пакет уже установлен на Вашем копьютере и загружен в рабочее пространство R; подробности вывода результатов анализа здесь не обсуждаются - см. справочные материалы по соответствующим командам).


# Пакет Hmisc, функция describe():
 
describe(mtcars)
mtcars 
 
 11  Variables      32  Observations
--------------------------------------------------------------------------------------------
mpg 
      n missing  unique    Mean     .05     .10     .25     .50     .75     .90     .95 
     31       1      25      20   11.85   14.30   15.35   19.20   22.15   30.40   31.40 
 
lowest : 10.4 13.3 14.3 14.7 15.0, highest: 26.0 27.3 30.4 32.4 33.9 
--------------------------------------------------------------------------------------------
cyl 
      n missing  unique    Mean 
     32       0       3   6.188 
 
4 (11, 34%), 6 (7, 22%), 8 (14, 44%) 
--------------------------------------------------------------------------------------------
disp 
      n missing  unique    Mean     .05     .10     .25     .50     .75     .90     .95 
     32       0      27   230.7   77.35   80.61  120.83  196.30  326.00  396.00  449.00 
 
lowest :  71.1  75.7  78.7  79.0  95.1, highest: 360.0 400.0 440.0 460.0 472.0 
--------------------------------------------------------------------------------------------
hp 
      n missing  unique    Mean     .05     .10     .25     .50     .75     .90     .95 
     32       0      22   146.7   63.65   66.00   96.50  123.00  180.00  243.50  253.55 
 
lowest :  52  62  65  66  91, highest: 215 230 245 264 335 
--------------------------------------------------------------------------------------------
drat 
      n missing  unique    Mean     .05     .10     .25     .50     .75     .90     .95 
     32       0      22   3.597   2.853   3.007   3.080   3.695   3.920   4.209   4.314 
 
lowest : 2.76 2.93 3.00 3.07 3.08, highest: 4.08 4.11 4.22 4.43 4.93 
--------------------------------------------------------------------------------------------
wt 
      n missing  unique    Mean     .05     .10     .25     .50     .75     .90     .95 
     32       0      29   3.217   1.736   1.956   2.581   3.325   3.610   4.048   5.293 
 
lowest : 1.513 1.615 1.835 1.935 2.140, highest: 3.845 4.070 5.250 5.345 5.424 
--------------------------------------------------------------------------------------------
qsec 
      n missing  unique    Mean     .05     .10     .25     .50     .75     .90     .95 
     32       0      30   17.85   15.05   15.53   16.89   17.71   18.90   19.99   20.10 
 
lowest : 14.50 14.60 15.41 15.50 15.84, highest: 19.90 20.00 20.01 20.22 22.90 
--------------------------------------------------------------------------------------------
vs 
      n missing  unique 
     32       0       2 
 
0 (18, 56%), 1 (14, 44%) 
--------------------------------------------------------------------------------------------
am 
      n missing  unique 
     32       0       2 
 
0 (19, 59%), 1 (13, 41%) 
--------------------------------------------------------------------------------------------
gear 
      n missing  unique    Mean 
     32       0       3   3.688 
 
3 (15, 47%), 4 (12, 38%), 5 (5, 16%) 
--------------------------------------------------------------------------------------------
carb 
      n missing  unique    Mean 
     32       0       6   2.812 
 
           1  2 3  4 6 8
Frequency  7 10 3 10 1 1
%         22 31 9 31 3 3
--------------------------------------------------------------------------------------------



# Пакет pastecs, функция stat.desc()
 
stat.desc(mtcars) 
                mpg    cyl    disp      hp    drat     wt   qsec
nbr.val       31.00  32.00 3.2e+01   32.00  32.000  32.00  32.00
nbr.null       0.00   0.00 0.0e+00    0.00   0.000   0.00   0.00
nbr.na         1.00   0.00 0.0e+00    0.00   0.000   0.00   0.00
min           10.40   4.00 7.1e+01   52.00   2.760   1.51  14.50
max           33.90   8.00 4.7e+02  335.00   4.930   5.42  22.90
range         23.50   4.00 4.0e+02  283.00   2.170   3.91   8.40
sum          620.10 198.00 7.4e+03 4694.00 115.090 102.95 571.16
median        19.20   6.00 2.0e+02  123.00   3.695   3.33  17.71
mean          20.00   6.19 2.3e+02  146.69   3.597   3.22  17.85
SE.mean        1.10   0.32 2.2e+01   12.12   0.095   0.17   0.32
CI.mean.0.95   2.24   0.64 4.5e+01   24.72   0.193   0.35   0.64
var           37.28   3.19 1.5e+04 4700.87   0.286   0.96   3.19
std.dev        6.11   1.79 1.2e+02   68.56   0.535   0.98   1.79
coef.var       0.31   0.29 5.4e-01    0.47   0.149   0.30   0.10
             vs am   gear  carb
nbr.val      NA NA  32.00 32.00
nbr.null     NA NA   0.00  0.00
nbr.na       NA NA   0.00  0.00
min          NA NA   3.00  1.00
max          NA NA   5.00  8.00
range        NA NA   2.00  7.00
sum          NA NA 118.00 90.00
median       NA NA   4.00  2.00
mean         NA NA   3.69  2.81
SE.mean      NA NA   0.13  0.29
CI.mean.0.95 NA NA   0.27  0.58
var          NA NA   0.54  2.61
std.dev      NA NA   0.74  1.62
coef.var     NA NA   0.20  0.57



# Пакет psych, функция describe()
 
describe(mtcars)
     var  n  mean     sd median trimmed    mad  min   max range
mpg    1 31  20.0   6.11   19.2    19.6   5.49 10.4  33.9  23.5
cyl    2 32   6.2   1.79    6.0     6.2   2.97  4.0   8.0   4.0
disp   3 32 230.7 123.94  196.3   222.5 140.48 71.1 472.0 400.9
hp     4 32 146.7  68.56  123.0   141.2  77.10 52.0 335.0 283.0
drat   5 32   3.6   0.53    3.7     3.6   0.70  2.8   4.9   2.2
wt     6 32   3.2   0.98    3.3     3.1   0.77  1.5   5.4   3.9
qsec   7 32  17.9   1.79   17.7    17.8   1.42 14.5  22.9   8.4
vs*    8 32   1.4   0.50    1.0     1.4   0.00  1.0   2.0   1.0
am*    9 32   1.4   0.50    1.0     1.4   0.00  1.0   2.0   1.0
gear  10 32   3.7   0.74    4.0     3.6   1.48  3.0   5.0   2.0
carb  11 32   2.8   1.62    2.0     2.6   1.48  1.0   8.0   7.0
      skew kurtosis    se
mpg   0.64    -0.03  1.10
cyl  -0.17    -1.76  0.32
disp  0.38    -1.07 21.91
hp    0.73     0.28 12.12
drat  0.27    -0.45  0.09
wt    0.42     0.42  0.17
qsec  0.37     0.86  0.32
vs*   0.24    -2.06  0.09
am*   0.36    -1.97  0.09
gear  0.53    -0.90  0.13
carb  1.05     2.02  0.29



# Пакет psych, функция describe.by() - расчет параметров описательной статистики
# для каждого уровня некоторого фактора:
 
describe.by(mtcars, mtcars$am)
group: 0
     var  n  mean     sd median trimmed    mad   min   max range
mpg    1 19  17.1   3.83   17.3    17.1   3.11  10.4  24.4  14.0
cyl    2 19   7.0   1.54    8.0     7.1   0.00   4.0   8.0   4.0
disp   3 19 290.4 110.17  275.8   289.7 124.83 120.1 472.0 351.9
hp     4 19 160.3  53.91  175.0   161.1  77.10  62.0 245.0 183.0
drat   5 19   3.3   0.39    3.1     3.3   0.22   2.8   3.9   1.2
wt     6 19   3.8   0.78    3.5     3.8   0.45   2.5   5.4   3.0
qsec   7 19  18.2   1.75   17.8    18.1   1.19  15.4  22.9   7.5
vs*    8 19   1.4   0.50    1.0     1.4   0.00   1.0   2.0   1.0
am*    9 19   1.0   0.00    1.0     1.0   0.00   1.0   1.0   0.0
gear  10 19   3.2   0.42    3.0     3.2   0.00   3.0   4.0   1.0
carb  11 19   2.7   1.15    3.0     2.8   1.48   1.0   4.0   3.0
      skew kurtosis    se
mpg   0.01    -0.33  0.88
cyl  -0.95    -0.24  0.35
disp  0.05    -1.01 25.28
hp   -0.01    -0.93 12.37
drat  0.50    -1.06  0.09
wt    0.98     1.06  0.18
qsec  0.85     1.66  0.40
vs*   0.50    -1.86  0.11
am*    NaN      NaN  0.00
gear  1.31     0.42  0.10
carb -0.14    -1.47  0.26
------------------------------------------------ 
group: 1
     var  n  mean    sd median trimmed   mad  min   max range
mpg    1 12  24.5  6.42   23.7    24.5  7.93 15.0  33.9  18.9
cyl    2 13   5.1  1.55    4.0     4.9  0.00  4.0   8.0   4.0
disp   3 13 143.5 87.20  120.3   131.2 58.86 71.1 351.0 279.9
hp     4 13 126.8 84.06  109.0   114.7 63.75 52.0 335.0 283.0
drat   5 13   4.0  0.36    4.1     4.0  0.27  3.5   4.9   1.4
wt     6 13   2.4  0.62    2.3     2.4  0.68  1.5   3.6   2.1
qsec   7 13  17.4  1.79   17.0    17.4  2.34 14.5  19.9   5.4
vs*    8 13   1.5  0.52    2.0     1.6  0.00  1.0   2.0   1.0
am*    9 13   2.0  0.00    2.0     2.0  0.00  2.0   2.0   0.0
gear  10 13   4.4  0.51    4.0     4.4  0.00  4.0   5.0   1.0
carb  11 13   2.9  2.18    2.0     2.6  1.48  1.0   8.0   7.0
      skew kurtosis    se
mpg  -0.01    -1.35  1.85
cyl   0.87    -0.15  0.43
disp  1.33     2.17 24.19
hp    1.36     2.46 23.31
drat  0.79     1.83  0.10
wt    0.21    -0.65  0.17
qsec -0.23    -1.10  0.50
vs*  -0.14    -2.36  0.14
am*    NaN      NaN  0.00
gear  0.42    -2.06  0.14
carb  0.98     1.07  0.60



# Пакет doBy. Обладает мощным функционалом. Пример приведен ниже:
 
summaryBy(mpg + wt ~ cyl + vs, data = mtcars,
          FUN = function(x) { c(m = mean(x), s = sd(x)) } )
  cyl vs mpg.m mpg.s wt.m wt.s
1   4  0    26    NA  2.1   NA
2   4  1    NA    NA  2.3 0.60
3   6  0    21  0.75  2.8 0.13
4   6  1    19  1.63  3.4 0.12
5   8  0    15  2.56  4.0 0.76



2 комментария :

Denis Atroshenko комментирует...

Здравствуйте. По-моему здесь небольшая опечатка в задании функции:
"Для начала создадим функцию SE для расчета стандартных ошибок (см. пример выше):

SE <- function(x) {sd(x)/length(sqrt(x))}

Теперь совместим эту новую функцию с tapply():"

Sergey Mastitsky комментирует...

Спасибо, Денис! Перепутал функции местами, но сейчас исправил.

Отправить комментарий