сабсетов из объекта phyloseq?
Функция для сабсетов для таких объектов основана на обычном subset() для датафреймов (код в самом низу). Я бы хотела применить эту функцию к списку выражений (который у меня формируется автоматически из названия фактора и нужных его уровней), но что-то идет не так...
Какой может быть самый простой обходной путь? Написать на основе их функции свою multisubset_samples(), которая возвращает список сабсетов?
# Обычный subset для data.frame работает
dfr <- data.frame(group = c("A", "A", "B", "A", "B"), value = 1:5)
conditions <- list(expression(group == "A" & value < 3),
expression(group == "A"),
expression(group == "B"))
lapply(conditions, function(x) subset(dfr, eval(x)))
# Создаем объект phyloseq
# if (!requireNamespace("BiocManager", quietly = TRUE))
# install.packages("BiocManager")
# BiocManager::install("phyloseq")
library("phyloseq")
OTU <- otu_table(matrix(1:20, ncol = 4), taxa_are_rows = TRUE)
samples <- sample_data(dfr)
mb <- phyloseq(OTU, samples)
# С одним условием все работает
subset_samples(mb, group == "A")
# Со списком условий - нет
lapply(conditions, function(x) subset_samples(mb, eval(x)))
# Error in eval(x) : object 'value' not found
# Исходный код `subset_samples`
# function (physeq, ...)
# {
# if (is.null(sample_data(physeq))) {
# cat("Nothing subset. No sample_data in physeq.\n")
# return(physeq)
# }
# else {
# oldDF <- as(sample_data(physeq), "data.frame")
# newDF <- subset(oldDF, ...)
# if (class(physeq) == "sample_data") {
# return(sample_data(newDF))
# }
# else {
# sample_data(physeq) <- sample_data(newDF)
# return(physeq)
# }
# }
# }
UPD:
Пока написала корявую функцию для сплита по фактору.
library("onehot")
multisubset_samples <- function (physeq, fact)
{
if (is.null(sample_data(physeq))) {
cat("Nothing subset. No sample_data in physeq.\n")
return(physeq)
}
else {
groups <- levels(fact)
dfr <- data.frame(fact)
fltr <- predict(onehot(dfr), dfr)
colnames(fltr) <- (gsub("fact=", "", colnames(fltr)))
fltr <- as.list(data.frame(fltr>0))
oldDF <- as(sample_data(physeq), "data.frame")
newDF_list <- lapply(fltr, function(x) subset(oldDF, eval(x)))
if (class(physeq) == "sample_data") {
res <- lapply(newDF_list, sample_data)
return(res)
}
else {
n <- length(fltr)
res <- vector(mode = "list", length = n)
for(i in 1:n) {
res[[i]] <- physeq
sample_data(res[[i]]) <- sample_data(newDF_list[[i]])
}
return(res)
}
}
}
А может лучше формировать не выражения, а сразу уже вектор с TRUE, FALSE? Тогда потом уже просто выбирать подвыборки самыми стандартными средствами
там s4 класс, а не стандартный датафрейм. пакетная функция сабсета его внутри преобразует и сплитит, но прямым путем оно не очень решается
Это я понимаю. В этом и затык...
проблема не в S4, а в subset
В целом, это пример кода как нельзя делать. Видно, что писали не разработчики. Марина, к Вам вопросов вроде как нет, выглядите как жертва. Раскладываю по полочкам. 1. subset использует NSE, причем достаточно криво, ибо subset рассчитан был исключительно на интерактивное применение. 2. для программного использования был придумал tidyselect & rlang 3. поскольку тут целая экосистема, написанная подобным образом, то тратить время на разбор идеологически кривого кода нет ни малейшего смысла. подойдет любой паллиатив, включая лом. 4. Лезть в движок NSE не понимаю его основ и глубин — категорически запрещено. Использовать отлаженные решения - ДА, мастерить самому — НЕТ. Вот такая грубая конструкция работает так, как изначально хотелось: library(phyloseq) library(dplyr) dfr <- data.frame(group = c("A", "A", "B", "A", "B"), value = 1:5) conditions <- list(quote(group == "A" & value < 3), quote(group == "A"), quote(group == "B")) OTU <- otu_table(matrix(1:20, ncol = 4), taxa_are_rows = TRUE) samples <- sample_data(dfr) mb <- phyloseq(OTU, samples) lapply(conditions, function(x){eval(parse(text = paste0("subset_samples(mb, ", deparse(x),")")))}) По-хорошему, за такой код в продуктиве отрывают руки или голову. Зависит от размера катастрофы. Кому интересно понять, почему subset является страшным косяком — можете почитать литературу: - Why base::subset is bad? - [Is R base::subset() really that bad?](https://win-vector.com/2018/02/23/is-r-basesubset-really-that-bad/) - [Why is [ better than subset?](https://stackoverflow.com/questions/9860090/why-is-better-than-subset). - [Non-standard evaluation](http://adv-r.had.co.nz/Computing-on-the-language.html). Читаем тут ссылку про проблемы subset в п. [Calling from another function]
дак subset только как резалка базы, а дальше все руками. в таком виде он прекрасно работает.
предлагаю прочесть все от начала до конца, особенно мнение мэтров. в нормальном коде subset категорически запрещен к использованию.
Хм, но Хэдли ведь не икона языка. Там полно других именитых персон.
я привел мнение нескольких человек John Mount тоже весьма известная личность. + описание механики `subset`и примеры косяков по ссылкам все разобрано по косточкам
Про косяки с subset знаю. Со временем начну побеждать subset и в чужом коде.
Зачем эти ужасные parse/deparse, когда результат quote - уже разобранное выражение? Вдобавок deparse тупо сломается, если выражение слишком длинное, больше 60 символов - разобьет его на вектор с несколькими элементами. Правильный вариант такой (отличие только в последней строке): library(phyloseq) library(dplyr) dfr <- data.frame(group = c("A", "A", "B", "A", "B"), value = 1:5) conditions <- list(quote(group == "A" & value < 3), quote(group == "A"), quote(group == "B")) OTU <- otu_table(matrix(1:20, ncol = 4), taxa_are_rows = TRUE) samples <- sample_data(dfr) mb <- phyloseq(OTU, samples) lapply(conditions, function(x) eval( substitute( subset_samples(mb, cond), list(cond = x) ) ) )
Полностью согласен. я нарисовал на скорую руку, исключительно для демонстрации. Ограничение у deparse есть, в данном случае все было вроде короткое. Вечер пятницы не самое лучшее время для таких глубин. Сказать, что силен в NSE и читаю без словаря — не могу. Проблематика, тянущаяся еще со времен LISP ясна, но при прочих равных предпочитаю писать простой и понятный код без персонального NSE. Но смысл этого чата и состоит в том, чтобы сообща находить лучшее решение. Кооперативка. За ответ спасибо, супер!
deparse хорошо контролируется входными параметрами. https://www.rdocumentation.org/packages/base/versions/3.6.2/topics/deparse width.cutoff управляет длиной выхода, так что 60 по умолчанию некритично. просто дополнение.
Обсуждают сегодня