Какие выборки считаются связанными |
Здравствуйте, гость ( Вход | Регистрация )
Какие выборки считаются связанными |
2.12.2011 - 07:39
Сообщение
#1
|
|
Группа: Пользователи Сообщений: 33 Регистрация: 31.07.2008 Пользователь №: 5185 |
Всем добрый день! Столкнулась с проблемой и просто зациклилась на ней. Пример из книги В.Ю. Урбаха : две делянки пшеницы, одна опыт, вторая контроль, измерялась урожайность раз в год.
Год| 1947|1948|1949|1950|1951|1952|1953 Опыт|22.9|20.2|19.5|30.5|35.6|31.9|27.7 Контроль|19.4|16.2|16.9|29.3|31.4|28.5|26.6 Для сравнения урожайности применяется критерий Стьюдента для парных выборок. У меня аналогичная задача, но я не могу доказать, почему эти выборки следует считать парными. На все мои объяснения, что опыт и контроль связаны годом, и что нельзя сравнивать урожайность первой делянки за 1947 год и урожайность второй делянки, например, за 1953г., а необходимо рассматривать именно пары, мне рассказывают про пациентов до и после лечения, и что там да, связанные, а здесь никакой связи нет. Может быть я не права? А если права, то, как объяснить так, чтобы не у кого не возникало никаких сомнений. Помогите, пожалуйста! Заранее большое спасибо. |
|
2.12.2011 - 09:14
Сообщение
#2
|
|
Группа: Пользователи Сообщений: 1091 Регистрация: 26.08.2010 Пользователь №: 22699 |
У меня аналогичная задача, но я не могу доказать, почему эти выборки следует считать парными. Вы забыли сказать в чем заключается ваша задача. |
|
2.12.2011 - 10:29
Сообщение
#3
|
|
Группа: Пользователи Сообщений: 33 Регистрация: 31.07.2008 Пользователь №: 5185 |
Прошу прощения. У меня два города и показатели числа коек на 1000 населения в них за 5 лет:
Год|2006|2007|2008|2009|2010 Город 1|8.5|8.8|7.9|7.5|8.1 Город 2|10.6|12.4|10.9|10.6|7.7 |
|
2.12.2011 - 14:01
Сообщение
#4
|
|
Группа: Пользователи Сообщений: 1091 Регистрация: 26.08.2010 Пользователь №: 22699 |
Прошу прощения. У меня два города и показатели числа коек на 1000 населения в них за 5 лет: Год|2006|2007|2008|2009|2010 Город 1|8.5|8.8|7.9|7.5|8.1 Город 2|10.6|12.4|10.9|10.6|7.7 1 по моему сравнение с вегетацией растений и зависимостью от года не очень корректно в данном случае 2 модель данных другая, там все в течении года появляется и убирается. в случае койкомест большая их часть переходит от года к году. (может надо каким то образом отображать движение койкомест (прибыло - убыло))? 3 ну поскольку "глобус выдан" попробуем посмотреть на данные Код # считываем "как есть" > read.table("data.txt",sep="|") V1 V2 V3 V4 V5 V6 1 Год 2006.0 2007.0 2008.0 2009.0 2010.0 2 Город 1 8.5 8.8 7.9 7.5 8.1 3 Город 2 10.6 12.4 10.9 10.6 7.7 # преобразуем в "широкий" датафрейм data<-as.data.frame(t(read.table("data.txt",sep="|")[,2:6])) names(data)<-read.table("data.txt",sep="|")[,1] > data Год Город 1 Город 2 V2 2006 8.5 10.6 V3 2007 8.8 12.4 V4 2008 7.9 10.9 V5 2009 7.5 10.6 V6 2010 8.1 7.7 # преобразуем в "длинный" датафрейм library(reshape) > data.long<-melt(data=data,id.vars="Год", measure.vars=c("Город 1","Город 2")) > data.long Год variable value 1 2006 Город 1 8.5 2 2007 Город 1 8.8 3 2008 Город 1 7.9 4 2009 Город 1 7.5 5 2010 Город 1 8.1 6 2006 Город 2 10.6 7 2007 Город 2 12.4 8 2008 Город 2 10.9 9 2009 Город 2 10.6 10 2010 Город 2 7.7 # строим модель смешанных эффектов library(lme4) # полную > model.ful <- lmer(value ~ 1 + (1|Год) + (1|variable), data=data.long) > model.ful Linear mixed model fit by REML Formula: value ~ 1 + (1 | Год) + (1 | variable) Data: data.long AIC BIC logLik deviance REMLdev 41.91 43.12 -16.96 36 33.91 Random effects: Groups Name Variance Std.Dev. Год (Intercept) 0.3120 0.55857 variable (Intercept) 2.3455 1.53150 Residual 1.2685 1.12628 Number of obs: 10, groups: Год, 5; variable, 2 Fixed effects: Estimate Std. Error t value (Intercept) 9.300 1.167 7.969 # только год > model.год <- lmer(value ~ 1 + (1|Год) , data=data.long) > model.год Linear mixed model fit by REML Formula: value ~ 1 + (1 | Год) Data: data.long AIC BIC logLik deviance REMLdev 43.27 44.17 -18.63 37.79 37.27 Random effects: Groups Name Variance Std.Dev. Год (Intercept) 0.0000 0.0000 Residual 2.8489 1.6879 Number of obs: 10, groups: Год, 5 Fixed effects: Estimate Std. Error t value (Intercept) 9.3000 0.5337 17.42 #только город > model.город <- lmer(value ~ 1 + (1|variable), data=data.long) > model.город Linear mixed model fit by REML Formula: value ~ 1 + (1 | variable) Data: data.long AIC BIC logLik deviance REMLdev 40.07 40.98 -17.04 36.12 34.07 Random effects: Groups Name Variance Std.Dev. variable (Intercept) 2.2831 1.5110 Residual 1.5805 1.2572 Number of obs: 10, groups: variable, 2 Fixed effects: Estimate Std. Error t value (Intercept) 9.30 1.14 8.158 # отображаем все графически с доверительными интервалами dotplot(ranef(model.ful, data=data.long, postVar = TRUE))$Год dotplot(ranef(model.ful, data=data.long, postVar = TRUE))$variable dotplot(ranef(model.год, data=data.long, postVar = TRUE)) dotplot(ranef(model.город, data=data.long, postVar = TRUE)) # ну и сами данные тоже отображаем графически print(dotplot(reorder(Год, value) ~ value, data.long, groups = variable, ylab = "Год")) судя по всему год не имеет достоверного вклада в модель (в качестве фактора дизайна эксперимента), отмечается слабая тенденция. Книга (на простом английском языке) с описание методик анализа лежит здесь http://lme4.r-forge.r-project.org/book/ отличия между моделями Код > anova(model.ful,model.город)
Data: data.long Models: model.город: value ~ 1 + (1 | variable) model.ful: value ~ 1 + (1 | Год) + (1 | variable) Df AIC BIC logLik Chisq Chi Df Pr(>Chisq) model.город 3 42.116 43.024 -18.058 model.ful 4 44.004 45.215 -18.002 0.1121 1 0.7378 > anova(model.ful,model.год) Data: data.long Models: model.год: value ~ 1 + (1 | Год) model.ful: value ~ 1 + (1 | Год) + (1 | variable) Df AIC BIC logLik Chisq Chi Df Pr(>Chisq) model.год 3 43.794 44.702 -18.897 model.ful 4 44.004 45.215 -18.002 1.7903 1 0.1809 > anova(model.год,model.город) Data: data.long Models: model.год: value ~ 1 + (1 | Год) model.город: value ~ 1 + (1 | variable) Df AIC BIC logLik Chisq Chi Df Pr(>Chisq) model.год 3 43.794 44.702 -18.897 model.город 3 42.116 43.024 -18.058 1.6782 0 < 2.2e-16 *** --- Signif. codes: 0 ?***? 0.001 ?**? 0.01 ?*? 0.05 ?.? 0.1 ? ? 1 > Сообщение отредактировал p2004r - 2.12.2011 - 14:31 |
|
2.12.2011 - 14:26
Сообщение
#5
|
|
Группа: Пользователи Сообщений: 1091 Регистрация: 26.08.2010 Пользователь №: 22699 |
Введем год в модель в качестве предиктора
Код > model <- lmer(value ~ Год + (1|variable), data=data.long) > model Linear mixed model fit by REML Formula: value ~ Год + (1 | variable) Data: data.long AIC BIC logLik deviance REMLdev 39.45 40.66 -15.72 32.28 31.45 Random effects: Groups Name Variance Std.Dev. variable (Intercept) 2.3721 1.5402 Residual 1.1342 1.0650 Number of obs: 10, groups: variable, 2 Fixed effects: Estimate Std. Error t value (Intercept) 983.1800 478.1898 2.056 Год -0.4850 0.2381 -2.037 Correlation of Fixed Effects: (Intr) Год -1.000 Предупреждения 1: 'abbreviate' использована с не-ASCII символами 2: 'abbreviate' использована с не-ASCII символами > anova(model.ful,model) Data: data.long Models: model.ful: value ~ 1 + (1 | Год) + (1 | variable) model: value ~ Год + (1 | variable) Df AIC BIC logLik Chisq Chi Df Pr(>Chisq) model.ful 4 44.004 45.215 -18.002 model 4 40.284 41.494 -16.142 3.7201 0 < 2.2e-16 *** --- Signif. codes: 0 ?***? 0.001 ?**? 0.01 ?*? 0.05 ?.? 0.1 ? ? 1 > dotplot(ranef(model, data=data.long, postVar = TRUE)) города достоверно отличаются, влияние переменной "год" существенно. ну и отобразим графически эти простые модели Код print(xyplot(value ~ Год | variable, data.long, aspect = "xy",
panel = function(...) { panel.xyplot(...) panel.abline(lm(data.long$value[as.numeric(data.long$variable)==packet.number()] ~ data.long$Год[as.numeric(data.long$variable)==packet.number()])) panel.abline(fixef(model), col.line = "red", lty = 1 ) panel.abline(fixef(model.ful), col.line = "green", lty = 2 ) } )) Сообщение отредактировал p2004r - 2.12.2011 - 16:45 |
|
2.12.2011 - 15:03
Сообщение
#6
|
|
Группа: Пользователи Сообщений: 902 Регистрация: 23.08.2010 Пользователь №: 22694 |
Всем добрый день! Столкнулась с проблемой и просто зациклилась на ней. Пример из книги В.Ю. Урбаха : две делянки пшеницы, одна опыт, вторая контроль, измерялась урожайность раз в год. Год| 1947|1948|1949|1950|1951|1952|1953 Опыт|22.9|20.2|19.5|30.5|35.6|31.9|27.7 Контроль|19.4|16.2|16.9|29.3|31.4|28.5|26.6 Для сравнения урожайности применяется критерий Стьюдента для парных выборок. Stefa, скажите, а уважаемый Урбах описал вероятностно-статистическую модель порождения данных? Какая гипотеза проверяется? Если не описал - выбросьте эту книгу не только из головы, но и из дома. |
|
2.12.2011 - 19:48
Сообщение
#7
|
|
Группа: Пользователи Сообщений: 1091 Регистрация: 26.08.2010 Пользователь №: 22699 |
Год| 1947|1948|1949|1950|1951|1952|1953 Опыт|22.9|20.2|19.5|30.5|35.6|31.9|27.7 Контроль|19.4|16.2|16.9|29.3|31.4|28.5|26.6 Ну а для сельскохозяйственного опыта все в порядке, для дизайна эксперимента как раз значим год. А участок ничего не вносит. Код > data.s.long<-melt(data=data.s,id.vars="Год", measure.vars=c("Опыт","Контроль")) > data.s.long Год variable value 1 1947 Опыт 22.9 2 1948 Опыт 20.2 3 1949 Опыт 19.5 4 1950 Опыт 30.5 5 1951 Опыт 35.6 6 1947 Контроль 19.4 7 1948 Контроль 16.2 8 1949 Контроль 16.9 9 1950 Контроль 29.3 10 1951 Контроль 31.4 > model.s.ful <- lmer(value ~ 1 + (1|Год) + (1|variable), data=data.s.long) > model.s.ful Linear mixed model fit by REML Formula: value ~ 1 + (1 | Год) + (1 | variable) Data: data.s.long AIC BIC logLik deviance REMLdev 56.32 57.54 -24.16 52.62 48.32 Random effects: Groups Name Variance Std.Dev. Год (Intercept) 49.693 7.04933 variable (Intercept) 4.654 2.15731 Residual 0.755 0.86891 Number of obs: 10, groups: Год, 5; variable, 2 Fixed effects: Estimate Std. Error t value (Intercept) 24.190 3.513 6.886 > model.s.год <- lmer(value ~ 1 + (1|Год) , data=data.s.long) > model.s.год Linear mixed model fit by REML Formula: value ~ 1 + (1 | Год) Data: data.s.long AIC BIC logLik deviance REMLdev 60.71 61.62 -27.36 58.8 54.71 Random effects: Groups Name Variance Std.Dev. Год (Intercept) 47.366 6.8823 Residual 5.409 2.3257 Number of obs: 10, groups: Год, 5 Fixed effects: Estimate Std. Error t value (Intercept) 24.190 3.164 7.645 > model.s.опыт <- lmer(value ~ 1 + (1|variable), data=data.s.long) > model.s.опыт Linear mixed model fit by REML Formula: value ~ 1 + (1 | variable) Data: data.s.long AIC BIC logLik deviance REMLdev 68.59 69.5 -31.3 65.94 62.59 Random effects: Groups Name Variance Std.Dev. variable (Intercept) 0.000 0.0000 Residual 47.512 6.8929 Number of obs: 10, groups: variable, 2 Fixed effects: Estimate Std. Error t value (Intercept) 24.19 2.18 11.1 |
|
3.12.2011 - 16:10
Сообщение
#8
|
|
Группа: Пользователи Сообщений: 1202 Регистрация: 13.01.2008 Из: Челябинск Пользователь №: 4704 |
>100$. Вообще-то учебник Урбаха весьма неплох для своего времени и особенно - места. В нём есть почти уникальные штуки, которые по крупицам рассыпаны в другой литературе, а кое-чего в других книгах я просто не нашёл, хотя по идее быть должно. Также он содержит больше доживших до наших дней методов, чем его ровесники от Плохинского. Хотя именно по зависимым выборкам я тоже не согласен. См. ниже.
>Stefa. На картинке я привёл цитату из соответствующего места Урбаха (стр. 155). На мой взгляд здесь есть путаница двух понятий: истинной завимости групп и преобладающего источника изменчивости. Настоящие зависимые выборки - опыты на одних и тех же животных и опыты на делянках с одних и тех же полей в разные годы. В первой же части цитаты описывается ситуация, которая не подразумевает зависимого характером выборок: то, что для каких-то биологических признаков именно год, а не поле привносит большую изменчивость, которую нужно грамотно учесть, ещё не делает выборки зависимыми. Просто получается что год - более существенный для урожайности фактор, чем поле. Если мы будем сравнивать урожайность с разных полей для одного ряда лет, а поля будут разные - некорректно считать и учитывать только ошибку разности в соответвующих парах лет, т.к. остаётся неучтённой изменчивость самих полей. А значит - это не есть зависимые выборки. Зависимые были бы тогда, когда мы бы поделили несколько полей и на одной части удобряли, на другой - нет. Тогда можно было считать, что разность между значениями на одном поле отражает именно результат удобрения и в обобщающей разности была бы только один источник изменчивости - разная "реакция" полей на удобрение. Я не особо разбирался что сделал р2004r, но вашу задачу можно решить: (1) двухфакторным дисперсионным анализом с единственным наблюдением на ячейку комплекса (параметрический подход) или (2) анализом Фридмана или Квейд (непараметрический подход). >p2004r. Если работать в классическом русле параметрикой, то это - двухфакторный дисперсионным анализом с единственным наблюдением на ячейку комплекса. Фактор "Город" - фиксированный, фактор "Год" - случайный. Взаимодействие факторов "Город х Год" неотделимо от ошибки анализа (т.к. в ячейках только по одной цифре) и само выступает в качестве статистики ошибки для главных эффектов. Ниже я прикрепил результаты этого анализа. Город значим, год - незначим. Год "вкладывает" в изменчивость койко-мест существенно меньше (лень считать компоненты дисперсии) и не является существенным. |
|
3.12.2011 - 17:59
Сообщение
#9
|
|
Группа: Пользователи Сообщений: 902 Регистрация: 23.08.2010 Пользователь №: 22694 |
> Nokh, коль скоро в посте #3 не содержится даже минимальных сведений о задаче, которую
Кроме того, мне не совсем понятно вот что: здесь мы имеем дело с (очевидно) независимыми выборками (2 города), а Фридман и Квейд работают со связанными выборками. Их-то на кой советуете? И еще. Вводить год в модель в качестве регрессора в данном контексте-очевидная глупость. Время - признак, не имеющий градаций. Вот рассмотреть в качестве предиктора величину инвестиций в здравоохранение (по годам для каждого города) - другое дело. Тогда, возможно, станет понятно, почему в городе Б больше коек на 1000 населения, нежели в граде А. Кстати, раз уж Урбах у Вас под рукой, что он пишет в смысле (возможного) ответа на мой вопрос из поста #6? P.S. Не люблю ломать копья в отсутствие топик-стартера. Сообщение отредактировал 100$ - 3.12.2011 - 18:19 |
|
3.12.2011 - 18:27
Сообщение
#10
|
|
Группа: Пользователи Сообщений: 1202 Регистрация: 13.01.2008 Из: Челябинск Пользователь №: 4704 |
Кроме того, мне не совсем понятно вот что: здесь мы имеем дело с (очевидно) независимыми выборками (2 города), а Фридман и Квейд работают со связанными выборками. Их-то на кой советуете? Совсем необязательно с зависимыми. Это - непараметрические аналоги двухфакторного дисперсионного анализа с рандомизированными блоками. Ссылки искать не хочется, но и википедии есть про пригодность для полноблочных экспериментальных планов. Хотя применяются и для зависимых выборок. Кстати, раз уж Урбах у Вас под рукой, что он пишет в смысле (возможного) ответа на мой вопрос из поста #6? А ничего не пишет... Но книгу не выброшу! : |
|
3.12.2011 - 19:17
Сообщение
#11
|
|
Группа: Пользователи Сообщений: 902 Регистрация: 23.08.2010 Пользователь №: 22694 |
Совсем необязательно с зависимыми. Это - непараметрические аналоги двухфакторного дисперсионного анализа с рандомизированными блоками. Ссылки искать не хочется, но и википедии есть про пригодность для полноблочных экспериментальных планов. Хотя применяются и для зависимых выборок. Надо же, как далеко продвинулась наука с тех пор, как я занимался ею в последний раз. Т.е. простенький критерий знаков (sign test) применяется для поиска различий в параметрах положения для двух зависимых выборок, а тест Фридмана, роспространяющий эту идею на случай большего количества выборок, каким-то непостижимым образом стал непараметрическим аналогом аж двухфакторного ANOVA? P.S. Оригинал статьи Квейд 1979 года в Тырнете выложен. Но я смотрел только на формулы, желая подкорректировать Кобзарёвы опечатки. Придется читать полностью. Сообщение отредактировал 100$ - 3.12.2011 - 19:38 |
|
3.12.2011 - 19:24
Сообщение
#12
|
|
Группа: Пользователи Сообщений: 902 Регистрация: 23.08.2010 Пользователь №: 22694 |
|
|
3.12.2011 - 19:38
Сообщение
#13
|
|
Группа: Пользователи Сообщений: 1202 Регистрация: 13.01.2008 Из: Челябинск Пользователь №: 4704 |
Надо же...И потом, сами посудите... Мне полемика не интересна. Я, в отличие от Вас, не бросаю нервных фраз в направлении авторов книг, типикстартеров и участников форума, оставляя вопрос без ответа. Я зашёл и ответил на тот вопрос топикстартера, который за 7 постов так и остался без ответа. О качестве выбора смайликов - тоже не вам судить, т.к. мои смайлики отражают мои эмоции, которые совсем необязательно должны казаться адекватными отдельным участникам форума, адекватность которых также может вызывать сомнения у других его участников. Сообщение отредактировал nokh - 3.12.2011 - 19:41 |
|
3.12.2011 - 20:10
Сообщение
#14
|
|
Группа: Пользователи Сообщений: 902 Регистрация: 23.08.2010 Пользователь №: 22694 |
Тогда предлагаю нам обоим посты 12 и 13 удалить собственноручно. Или пусть модераторы расстараются. Обещаю больше не экспериментировать с вашим самолюбием.
Сообщение отредактировал 100$ - 3.12.2011 - 20:20 |
|
4.12.2011 - 16:07
Сообщение
#15
|
|
Группа: Пользователи Сообщений: 1091 Регистрация: 26.08.2010 Пользователь №: 22699 |
И еще. Вводить год в модель в качестве регрессора в данном контексте-очевидная глупость. Время - признак, не имеющий градаций. Вот рассмотреть в качестве предиктора величину инвестиций в здравоохранение (по годам для каждого города) - другое дело. Тогда, возможно, станет понятно, почему в городе Б больше коек на 1000 населения, нежели в граде А. Ну ну... "глупость" это лонгитудинальное наблюдение считать независимым от времени. И _никакой_ учтенный фактор "по годам" природу модели в части зависимости значения текущего года от предыдущего не изменит (да и сам скорее всего будет зависим от своих предыдущих по времени значений). |
|