Я прочитал много сообщений, связанных с обработкой данных и «повторным» t-тестом, но я не могу понять, как это сделать в моем случае.
Вы можете получить мой пример набора данных для StackOverflow здесь: //www.dropbox.com/s/0b618fs1jjnuzbg/dataset.example.stckovflw.txt?dl=0
У меня есть большой кадр данных выражения gen, например:
> b<-read.delim("dataset.example.stckovflw.txt")
> head(b)
animal gen condition tissue LogFC
1 animalcontrol1 kjhss1 control brain 7.129283
2 animalcontrol1 sdth2 control brain 7.179909
3 animalcontrol1 sgdhstjh20 control brain 9.353147
4 animalcontrol1 jdygfjgdkydg21 control brain 6.459432
5 animalcontrol1 shfjdfyjydg22 control brain 9.372865
6 animalcontrol1 jdyjkdg23 control brain 9.541097
> str(b)
'data.frame': 21507 obs. of 5 variables:
$ animal : Factor w/ 25 levels "animalcontrol1",..: 1 1 1 1 1 1 1 1 1 1 ...
$ gen : Factor w/ 1131 levels "dghwg1041","dghwg1086",..: 480 761 787 360 863 385 133 888 563 738 ...
$ condition: Factor w/ 5 levels "control","treatmentA",..: 1 1 1 1 1 1 1 1 1 1 ...
$ tissue : Factor w/ 2 levels "brain","heart": 1 1 1 1 1 1 1 1 1 1 ...
$ LogFC : num 7.13 7.18 9.35 6.46 9.37 ...
Каждая группа состоит из 5 животных, и у каждого животного количественно определено много родов. (Однако у каждого животного может быть различный набор количественных родов, но также многие роды будут общими между животными и группами).
Я хотел бы выполнить t-тест для каждого поколения между моей обработанной группой (A, B, C или D) и контролем. Данные должны быть представлены в виде таблицы, содержащей p-значение для каждого поколения в каждой группе.
Поскольку у меня так много родов (тысячи), я не могу подмножить каждый род.
Вы знаете, как я могу автоматизировать процедуру?
Я думал о цикле, но я абсолютно не уверен, что он может достичь того, чего я хочу, и как действовать дальше.
Кроме того, я больше просматривал эти сообщения, используя функцию apply: Применить t-тест ко многим столбцам в фрейме данных, разделенному по факторам и Перебор t.tests для подмножеств фреймов данных в r
@andrew_reece: Большое спасибо за это. Это почти то, что я искал. Однако я не могу найти способ сделать это с помощью t-теста. ANOVA — это интересная информация, но тогда мне нужно будет знать, какая из обработанных групп значительно отличается от моей контрольной группы. Также мне нужно было бы знать, какая обработанная группа значительно отличается друг от друга, «два на два».
Я пытался использовать ваш код, изменив «aov (..)» в «t.test (…)». Для этого сначала я реализую подмножество (b, условие == «контроль» | условие == «лечениеA»), чтобы сравнить только две группы. Однако при экспорте таблицы результатов в файл csv таблица непонятна (без имени гена, без p-значений и т. д., только числа). Я буду продолжать искать способ сделать это правильно, но пока я застрял.
@42:
Большое спасибо за эти советы. Это просто пример набора данных, давайте предположим, что нам нужно использовать отдельные t-тесты.
Это очень полезное начало для изучения моих данных. Например, я пытался представить свои данные с помощью Venndiagrams. Я могу написать свой код, но он немного выходит за рамки первоначальной темы. Кроме того, я не знаю, как менее придирчиво обобщить общий «ген», обнаруженный в каждой комбинации условий, поэтому я упростил только 3 условия.
# Visualisation of shared genes by VennDiagrams :
# let's simplify and consider only 3 conditions :
b<-read.delim("dataset.example.stckovflw.txt")
b<- subset(b, condition == "control" | condition == "treatmentA" | condition == "treatmentB")
b1<-table(b$gen, b$condition)
b1
b2<-subset(data.frame(b1, "control" > 2
|"treatmentA" > 2
|"treatmentB" > 2 ))
b3<-subset(b2, Freq>2) # select only genes that have been quantified in at least 2 animals per group
b3
b4 = within(b3, {
Freq = ifelse(Freq > 1, 1, 0)
}) # for those observations, we consider the gene has been detected so we change the value 0 regardless the freq of occurence (>2)
b4
b5<-table(b4$Var1, b4$Var2)
write.csv(b5, file = "b5.csv")
# make an intermediate file .txt (just add manually the name of the cfirst column title)
# so now we have info
bb5<-read.delim("bb5.txt")
nrow(subset(bb5, control == 1))
nrow(subset(bb5, treatmentA == 1))
nrow(subset(bb5, treatmentB == 1))
nrow(subset(bb5, control == 1 & treatmentA == 1))
nrow(subset(bb5, control == 1 & treatmentB == 1))
nrow(subset(bb5, treatmentA == 1 & treatmentB == 1))
nrow(subset(bb5, control == 1 & treatmentA == 1 & treatmentB == 1))
library(grid)
library(futile.logger)
library(VennDiagram)
venn.plot <- draw.triple.venn(area1 = 1005,
area2 = 927,
area3 = 943,
n12 = 843,
n23 = 861,
n13 = 866,
n123 = 794,
category = c("controls", "treatmentA", "treatmentB"),
fill = c("red", "yellow", "blue"),
cex = 2,
cat.cex = 2,
lwd = 6,
lty = 'dashed',
fontface = "bold",
fontfamily = "sans",
cat.fontface = "bold",
cat.default.pos = "outer",
cat.pos = c(-27, 27, 135),
cat.dist = c(0.055, 0.055, 0.085),
cat.fontfamily = "sans",
rotation = 1);
genиgensдействительно требует рассмотрения этих значений как гена, а гены как генетические. - person IRTFM schedule 13.05.2018@ссылки в вашем фактическом сообщении не помечаются в учетных записях пользователей. Я бы не увидел ваш обновленный запрос, если бы не проверял. Лучше публиковать ответы, относящиеся к конкретному ответу, в виде комментариев к данному ответу. - person andrew_reece schedule 13.05.2018