Пункты наименьших квадратов для групп факторов уровня

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

Пример данных:

set.seed(1) 
model.data <- data.frame(time = factor(paste0("day", rep(1:8, each = 16))), 
    animal = factor(rep(1:16, each = 8)), 
    tissue = factor(c("blood", "liver", "kidney", "brain")), 
    value = runif(128) 
)

Настройка пользовательских контрасты для фактора "времени":

library("phia") 
custom.contrasts <- as.data.frame(contrastCoefficients(
    time ~ (day1+day2+day3)/3 - (day4+day5+day6)/3, 
    time ~ (day1+day2+day3)/3 - (day7+day8)/2, 
    time ~ (day4+day5+day6)/3 - (day7+day8)/2, 
    data = model.data, normalize = FALSE)) 

colnames(custom.contrasts) <- c("early - late", 
    "early - very late", 
    "late - very late") 

custom.contrasts.lsmc <- function(...) return(custom.contrasts)

Установка модели и вычисление наименьших квадратов означает:

library("lme4") 
tissue.model <- lmer(value ~ time * tissue + (1|animal), model.data) 
library("lsmeans") 
tissue.lsm <- lsmeans(tissue.model, custom.contrasts ~ time | tissue)

черчения :

plot(tissue.lsm$lsmeans) 
dev.new() 
plot(tissue.lsm$contrasts)

Теперь у второго участка есть комбинации, которые я хочу, но это показывает разницу между комбинированными средствами, а не средствами.

Я могу получить индивидуальные значения от tissue.lsm$lsmeans и рассчитать комбинированные средства самостоятельно, но у меня есть ощущение нытье, что есть более простой способ, которого я просто не вижу. В конце концов, все данные должны быть в lsmobj.

early.mean.liver = mean(model.data$value[model.data$tissue == "liver" & 
    model.data$time %in% c("day1", "day2", "day3")]) 
late.mean.liver = mean(model.data$value[model.data$tissue == "liver" & 
    model.data$time %in% c("day4", "day5", "day6")]) 
vlate.mean.liver = mean(model.data$value[model.data$tissue == "liver" & 
    model.data$time %in% c("day7", "day8")]) 
# ... for each level of "tissue" 


#compare to tissue.lsm$contrasts 
early.mean.liver - late.mean.liver 
early.mean.liver - vlate.mean.liver 
late.mean.liver - vlate.mean.liver

Я с нетерпением жду ваших замечаний или предложений. Благодаря!

источник

2015-08-02 tractorjack

Может просто создать второй коэффициенты матрицы, чтобы получить комбинированные средства, представляющие интерес для каждой ткани? – aosmith

Благодарим вас за предложение @aosmith! Я не знаю, как это сделать, хотя, можете ли вы привести мне пример, пожалуйста? – tractorjack

Один из вариантов заключается в вычислении коэффициентов контрастности для интересующего группового средства в дополнение к коэффициентам контрастности для различий в групповом средстве, который вы рассчитали в custom_contrasts. Например, вы можете сделать это отдельно как custom.contrasts2.

custom.contrasts2 <- as.data.frame(contrastCoefficients(
    time ~ (day1+day2+day3)/3, 
    time ~ (day4+day5+day6)/3, 
    time ~ (day7+day8)/2, 
    data = model.data, normalize = FALSE)) 

colnames(custom.contrasts2) <- c("early", 
          "late", 
          "very late") 

custom.contrasts2.lsmc <- function(...) return(custom.contrasts2) 

lsmeans(tissue.model, custom.contrasts2 ~ time | tissue)$contrasts

Вот только выходы для liver, которые группа означает, что вы после этого.

... 
tissue = liver: 
contrast estimate   SE df t.ratio p.value 
early  0.4481244 0.07902715 70.4 5.671 <.0001 
late  0.4618041 0.07902715 70.4 5.844 <.0001 
lvery late 0.3824247 0.09678810 70.4 3.951 0.0002

Если вы знаете, вы хотите, как группа означает и различие в группе означает, что вы можете просто добавить контраст коэффициентов матрицы создаются с помощью contrastCoefficients.

custom.contrasts <- as.data.frame(contrastCoefficients(
    time ~ (day1+day2+day3)/3, 
    time ~ (day4+day5+day6)/3, 
    time ~ (day7+day8)/2, 
    time ~ (day1+day2+day3)/3 - (day4+day5+day6)/3, 
    time ~ (day1+day2+day3)/3 - (day7+day8)/2, 
    time ~ (day4+day5+day6)/3 - (day7+day8)/2, 
    data = model.data, normalize = FALSE))

А затем назовите и выполните функцию .lsmc.

источник

2015-08-04 20:57:23 aosmith

Это именно то, что я искал! Спасибо огромное! Я понятия не имел, что могу использовать такие контрасты. – tractorjack

Следуя примеру @ aosmith в:

custom.means <- as.data.frame(contrastCoefficients(
    time ~ (day1+day2+day3)/3, 
    time ~ (day4+day5+day6)/3, 
    time ~ (day7+day8)/2, 
    data = model.data, normalize = FALSE)) 

colnames(custom.means) <- c("early", 
    "late", 
    "very late") 

custom.means.lsmc <- function(...) return(custom.means) 

tissue.means <- confint(lsmeans(tissue.model, custom.means ~ time | tissue)$contrasts) 

library("ggplot2") 
p <- ggplot(tissue.means, 
    aes(x = contrast, y = estimate, ymin = lower.CL, ymax = upper.CL)) + 
    geom_errorbar() + facet_wrap(~ tissue, ncol = 4) + xlab("time") 

print(p)

источник

2015-08-04 23:53:09 tractorjack

Пункты наименьших квадратов для групп факторов уровня

ответ

Смежные вопросы