Петля через t.тесты для подмножеств фреймов данных в r

У меня есть фрейм данных ‘math.numeric ‘ с 32 переменными. Каждая строка представляет учащегося, а каждая переменная является атрибутом. Студенты были разделены на 5 групп в зависимости от их итоговой оценки.

Данные выглядят следующим образом:

head(math.numeric)
school sex age address famsize Pstatus Medu Fedu Mjob Fjob reason ... group
1      1   18  2       1       1       4    4    1    5    1          2
1      1   17  2       1       2       1    1    1    3    1          2
1      1   15  2       2       2       1    1    1    3    3          3
1      1   15  2       1       2       4    2    2    4    2          4
1      1   16  2       1       2       3    3    3    3    2          3
1      2   16  2       2       2       4    3    4    3    4          4

Я выполняю t-тесты по каждой переменной для группы 1 против всех других групп, чтобы идентифицировать значительно отличающиеся атрибуты с этой группой. Я ищу, чтобы вытащить p-значения для каждого теста, такие как:

t.test(subset(math.numeric$school, math.numeric$group == 1),
      subset(math.numeric$school, math.numeric$group != 1))$p.value
t.test(subset(math.numeric$sex, math.numeric$group == 1), 
        subset(math.numeric$sex, math.numeric$group != 1))$p.value
t.test(subset(math.numeric$age, math.numeric$group == 1), 
        subset(math.numeric$age, math.numeric$group != 1))$p.value

Я пытаюсь понять, как я могу создать цикл, чтобы сделать это вместо того, чтобы выписывать каждый тест по одному за раз. Я попробовал петлю for и lapply, но до сих пор мне не повезло.

Я довольно новичок в этом, поэтому любая помощь была бы оценена.

Кортни

3 ответа

  1. Как это?

    pvals <- numeric() #the vector of p values
    k <- 1 #in case you choose to use a subset not continuing from 1
    
    # "for(i in seq(1,5))" is just doing the pvalues for the first 5 columns. You could do a 
    # list, like "c(1,2,4)" (in place of "seq(1,5)"), to do tests for columns 1, 2, and 4. 
    # To do all of the columns, try "for(i in seq(1,(ncol(math.numeric)-1)))".
    
    for(i in seq(1,5)){
    
      # using your code to grab the p-values and store them in the kth element of "pvals"
      pvals[k] <- t.test(subset(math.numeric[,i], math.numeric$group == 1),
             subset(math.numeric[,i], math.numeric$group != 1))$p.value    
    
      #iterating the "pvals" vector entry counter
      k=k+1
    }
    pvals #printing the p values for each test
    
  2. Данных вашего примера недостаточно для проведения t-тестов по всем подгруппам. По этой причине я беру irisнабор данных, который содержит 3 вида растений: Setosa, Versicolor и Virginica. Это мои группы. Вам придется соответствующим образом настроить код. Ниже я покажу, как протестировать одну группу против всех других групп, одну группу против друг друга и все комбинации отдельных групп.

    Одна группа Против всех других групп вместе взятых:

    Во-первых, предположим, что я хочу сравнить Versicolor и Virginica с Setosa, т. е. Setosa-это мойgroup 1, с которым должны сравниваться все другие группы. Простой способ достичь желаемого-это:

    sapply(names(iris)[-ncol(iris)], function(x){
                 t.test(iris[iris$Species=="setosa", x], 
                        iris[iris$Species!="setosa", x])$p.value 
                        })
    Sepal.Length  Sepal.Width Petal.Length  Petal.Width 
    7.709331e-32 1.035396e-13 1.746188e-69 1.347804e-60 
    

    Здесь я задал имена различных переменных в наборе names(iris)данных, включая столбец, указывающий на группирующую переменную [-ncol(iris)](поскольку это последний столбец), в качестве вектора tosapply, который передает соответствующие имена в качестве аргументов функции, которую я определил.

    Одна группа против каждой из других групп:

    В случае, если вы хотите сделать групповые сравнения для всех групп, может быть полезно следующее: Во-первых, создайте фрейм данных всех комбинаций переменных группы x, которые вы собираетесь сделать, исключая, конечно, саму группирующую переменную и справочную группу. Это может быть достигнуто путем:

    comps <- expand.grid(unique(iris$Species)[-1], # excluding Setosa as reference group
                         names(iris)[-ncol(iris)] # excluding group column
                         )
    head(comps)
            Var1         Var2
    1 versicolor Sepal.Length
    2  virginica Sepal.Length
    3 versicolor  Sepal.Width
    4  virginica  Sepal.Width
    5 versicolor Petal.Length
    6  virginica Petal.Length
    

    Здесь Var1представлены различные виды и Var2различные переменные, для которых должны проводиться сравнения. Ссылка group 1или Setosa неявно в этом случае. Теперь я могу использовать applyдля создания тестов. Я делаю это, используя каждую строку compsаргумента as с двумя элементами, первый из которых указывает, какой поворот группы это, а второй аргумент указывает, какую переменную следует сравнить. Они будут использоваться для подмножества исходного фрейма данных.

    comps$pval <- apply(comps, 1, function(x) {
        t.test(iris[iris$Species=="setosa", x[2]], iris[iris$Species==x[1], x[2]])$p.value 
        } )
    

    где группа 1 aka Setosa жестко закодирована в функции. Это дает мне фрейм данных с p-значениями для всех комбинаций (с Setosa в качестве контрольной группы), так что они легко искать:

    head(comps)
            Var1         Var2         pval
    1 versicolor Sepal.Length 3.746743e-17
    2  virginica Sepal.Length 3.966867e-25
    3 versicolor  Sepal.Width 2.484228e-15
    4  virginica  Sepal.Width 4.570771e-09
    5 versicolor Petal.Length 9.934433e-46
    6  virginica Petal.Length 9.269628e-50
    

    Все комбинации групп:

    Вы можете легко развернуть вышеизложенное, чтобы создать фрейм данных, содержащий p-значения t-тестов для каждой комбинации групп. Один подход был бы:

    comps <- expand.grid(unique(iris$Species), unique(iris$Species), names(iris)[-ncol(iris)])
    

    Это теперь имеет три колонки. Первые две группы, а третья переменная:

    head(comps)
            Var1       Var2         Var3
    1     setosa     setosa Sepal.Length
    2 versicolor     setosa Sepal.Length
    3  virginica     setosa Sepal.Length
    4     setosa versicolor Sepal.Length
    5 versicolor versicolor Sepal.Length
    6  virginica versicolor Sepal.Length
    

    Вы можете использовать это для выполнения тестов:

    comps$pval <- apply(comps, 1, function(x) {
      t.test(iris[iris$Species==x[1], x[3]], iris[iris$Species==x[2], x[3]])$p.value 
    } )
    

    Я получаю сообщение об ошибке: что я должен сделать?

    t.test может выдавать сообщение об ошибке, если размер выборки слишком мал или значения постоянны для одной группы. Это проблематично, так как это может произойти только для определенных групп, и вы можете не знать заранее, какой из них. Тем не менее ошибка нарушит весь вызов функции apply, и вы не сможете увидеть какие-либо результаты.

    Способ обойти это и определить проблемные группы-обернуть функцию t.testdplyr::failwith(см. Также ?tryCatch). Чтобы показать, как это работает, рассмотрим следующее:

    smalln <- data.frame(a=1, b=2)
    t.test(smalln$a, smalln$b)
    > Error in t.test.default(smalln$a, smalln$b) : not enough 'x' observations
    
    failproof.t <- failwith(default="Some default of your liking", t.test, quiet = T)
    failproof.t(smalln$a, smalln$b)
    [1] "Some default of your liking"
    

    Таким образом, всякий t.testраз, когда будет выбрасывать ошибку, вы получаете символ в результате вместо этого, и вычисление продолжается с другими группами. Излишне говорить, что вы также можете установить defaultчисло или что-нибудь еще. Это не обязательно должен быть персонаж.

    Статистический отказ от ответственности:
    Сказав все это, заметьте, что проведение нескольких t-тестов не обязательно является хорошей статистической практикой. Вы можете настроить p-значения для учета нескольких тестов или использовать альтернативные процедуры тестирования, которые проводят совместные тесты.

  3. Рассмотрите возможность разделения фрейма данных по группам и использования mapply()по столбцам. Выходные данные становятся скомпилированным списком статистики тестов: Statistica, parameter, p-value, confid. интервал, etc.

    # FILTER ROWS AND SUBSET NUMERIC COLS
    group1df <- df[df$group==1, 1:ncol(df)-1]
    othgroupdf <- df[df$group!=1, 1:ncol(df)-1]
    
    # T-TEST FCT
    tfct <- function(v1, v2){
          t.test(v1, v2) 
    }
    
    # RUN T-TESTS BY COL, SAVE RESULTS TO LIST
    ttests <- mapply(tfct, group1df, othgroupdf)