Помогите, пожалуйста, написать алгоритм для R |
Здравствуйте, гость ( Вход | Регистрация )
Помогите, пожалуйста, написать алгоритм для R |
29.08.2012 - 12:02
Сообщение
#1
|
|
Группа: Пользователи Сообщений: 43 Регистрация: 4.01.2012 Пользователь №: 23400 |
Добрый день.
Помогите, пожалуйста, написать алгоритм для R. У меня две базы данных, одна из них включает номер пациента, дату каждого измерения (как правило, измерение проводилось в среднем каждуjу неделю в течении 1 года ) и номер недели (1,2,3...54), вторая база включает информацию о ежедневной температуре воздуха. Задача: нужно посчитать среднюю температуру, предшествующую каждому измерению. Например, для первой недели нужно рассчитать среднюю темп за 7 дней до измерения (28apr1999-7дней), для второй недели за период 28apr1999-07may1999 и т.д. FID studydate week 1 28apr1999 1 1 07may1999 2 1 14may1999 3 1 21may1999 4 Second dataset temp date 5.37 1999-04-20 2.13 1999-04-21 1.6 1999-04-22 -0.17 1999-04-23 2.53 1999-04-24 Спасибо! |
|
29.08.2012 - 16:37
Сообщение
#2
|
|
Группа: Пользователи Сообщений: 1091 Регистрация: 26.08.2010 Пользователь №: 22699 |
Добрый день. Помогите, пожалуйста, написать алгоритм для R. У меня две базы данных, одна из них включает номер пациента, дату каждого измерения (как правило, измерение проводилось в среднем каждуjу неделю в течении 1 года ) и номер недели (1,2,3...54), вторая база включает информацию о ежедневной температуре воздуха. Задача: нужно посчитать среднюю температуру, предшествующую каждому измерению. Например, для первой недели нужно рассчитать среднюю темп за 7 дней до измерения (28apr1999-7дней), для второй недели за период 28apr1999-07may1999 и т.д. FID studydate week 1 28apr1999 1 1 07may1999 2 1 14may1999 3 1 21may1999 4 Second dataset temp date 5.37 1999-04-20 2.13 1999-04-21 1.6 1999-04-22 -0.17 1999-04-23 2.53 1999-04-24 Спасибо! Пожалуйста! Я предлагаю посчитать скользящим окном среднее для всех дней в базе. При импорте формата даты первого файла на некоторых системах возможно понадобится переключить локаль (мне в utf8 понадобилось). Примеры слишком короткие для демонстрации работоспособности. Код > temp.time<-read.table("time.txt", header=TRUE)
> temp.time$date<-as.Date(temp.time$date) > str(temp.time) 'data.frame': 5 obs. of 2 variables: $ temp: num 5.37 2.13 1.6 -0.17 2.53 $ date:Class 'Date' num [1:5] 10701 10702 10703 10704 10705 > > library(caTools) > temp.time <- data.frame(temp.time, avgtemp=runmean(temp.time$temp, 3, alg= "exact", endrule="mean", align = "right")) > temp.time temp date avgtemp 1 5.37 1999-04-20 3.750000 2 2.13 1999-04-21 3.750000 3 1.60 1999-04-22 3.033333 4 -0.17 1999-04-23 1.186667 5 2.53 1999-04-24 1.320000 > lct <- Sys.getlocale("LC_TIME"); Sys.setlocale("LC_TIME", "C") > data<-read.table("data.txt", header=TRUE) > data$studydate<-as.Date(data$studydate, "%d%b%Y") > data FID studydate week 1 1 1999-04-28 1 2 1 1999-05-07 2 3 1 1999-05-14 3 4 1 1999-05-21 4 > str(data) 'data.frame': 4 obs. of 3 variables: $ FID : int 1 1 1 1 $ studydate:Class 'Date' num [1:4] 10709 10718 10725 10732 $ week : int 1 2 3 4 > Sys.setlocale("LC_TIME", lct) > temp.time$avgtemp[temp.time$date==data[1,2]] |
|
29.08.2012 - 18:43
Сообщение
#3
|
|
Группа: Пользователи Сообщений: 43 Регистрация: 4.01.2012 Пользователь №: 23400 |
Спасибо!
У меня вопрос, почему к(width of moving window)=3? Не совсем поняла, как работает предложенный вами метод.. Я думала, что нужно написать loop для решения этой задачи. Но мои попытки не дали желаемого результата.. airpollution <- read.table("C:.../Pay.csv",blank.lines.skip = TRUE,fill=T,header=T,sep=";") summary(airpollution) dat <- paste(airpollution[,3],"-",airpollution[,2],"-",airpollution[,1],sep="") dat <- as.Date(dat, format="%Y-%m-%d") pay <- cbind(airpollution,dat) x <- read.dta("C:/.../symptoms.dta") x$stday <- as.Date(x$stday, origin="1960-01-01",format="%Y-%m-%d") temp<- rep(NA, length(x$FID)) index <- 1:length(pay$dat) pm10<- rep(NA, length(x$FID)) for (i in 1:length(x$FID)) { start <- x$stday[i]==pay$dat s <- index[start] end <- x$stday[i]+1==pay$dat e <- index[end] pm10[i] <- mean(pay$pm10[s:e],na.rm = T) paytemp[i] <- mean(pay$temp[s:e],na.rm = T) } Не уверена насчет этой части end <- x$stday[i]+1==pay$dat +1 это значит, что следующая строка, или я ошибаюсь? |
|
29.08.2012 - 20:57
Сообщение
#4
|
|
Группа: Пользователи Сообщений: 1091 Регистрация: 26.08.2010 Пользователь №: 22699 |
Спасибо! У меня вопрос, почему к(width of moving window)=3? Не совсем поняла, как работает предложенный вами метод.. Я думала, что нужно написать loop для решения этой задачи. Но мои попытки не дали желаемого результата.. airpollution <- read.table("C:.../Pay.csv",blank.lines.skip = TRUE,fill=T,header=T,sep=";") summary(airpollution) dat <- paste(airpollution[,3],"-",airpollution[,2],"-",airpollution[,1],sep="") dat <- as.Date(dat, format="%Y-%m-%d") pay <- cbind(airpollution,dat) x <- read.dta("C:/.../symptoms.dta") x$stday <- as.Date(x$stday, origin="1960-01-01",format="%Y-%m-%d") temp<- rep(NA, length(x$FID)) index <- 1:length(pay$dat) pm10<- rep(NA, length(x$FID)) for (i in 1:length(x$FID)) { start <- x$stday[i]==pay$dat s <- index[start] end <- x$stday[i]+1==pay$dat e <- index[end] pm10[i] <- mean(pay$pm10[s:e],na.rm = T) paytemp[i] <- mean(pay$temp[s:e],na.rm = T) } Не уверена насчет этой части end <- x$stday[i]+1==pay$dat +1 это значит, что следующая строка, или я ошибаюсь? 1) на таких коротких данных я выбрал усреднение по 3 дням (или меньше за счет выбранного алгоритма) чтобы хоть что то работало для отладки. Программа работает просто --- для всех дней справочника температуры рассчитывается среднее за предыдущие 7 дней (у меня для отладки на Ваших лапидарных примерах таблиц -- за 3 дня). После этого Вы просто запрашиваете последней строкой из справочника среднюю температуру для даты инцидента. Скажите куда Вам это среднее надо подставить? И приведите более длинные куски данных. Иначе мне непонятно чем Вас не устраивает предложенное мной решение. 2) Ну что могу сказать, как говорится --- " Хороший программист на фортране, на любом языке напишет программу на фортране" . По поводу указанной строчки кода -- это какая то невозможная конструкция Я не использую циклов в R, и Вам не советую. В R все делается через *apply() операторы. Для использование for() нужна крайне исключительная причина. Прочитайте пожалуйста руководство "Введение в R" (или на английском, или на русском со страницы перевода http://m7876.wiki.zoho.com/Introduction-to-R.html , а то у меня какое то чувство что я его зря переводил) |
|
29.08.2012 - 21:17
Сообщение
#5
|
|
Группа: Пользователи Сообщений: 43 Регистрация: 4.01.2012 Пользователь №: 23400 |
1) на таких коротких данных я выбрал усреднение по 3 дням (или меньше за счет выбранного алгоритма) чтобы хоть что то работало для отладки. Программа работает просто --- для всех дней справочника температуры рассчитывается среднее за предыдущие 7 дней (у меня для отладки на Ваших лапидарных примерах таблиц -- за 3 дня). После этого Вы просто запрашиваете последней строкой из справочника среднюю температуру для даты инцидента. Скажите куда Вам это среднее надо подставить? И приведите более длинные куски данных. Иначе мне непонятно чем Вас не устраивает предложенное мной решение. 2) Ну что могу сказать, как говорится --- " Хороший программист на фортране, на любом языке напишет программу на фортране" . По поводу указанной строчки кода -- это какая то невозможная конструкция Я не использую циклов в R, и Вам не советую. В R все делается через *apply() операторы. Для использование for() нужна крайне исключительная причина. Прочитайте пожалуйста руководство "Введение в R" (или на английском, или на русском со страницы перевода http://m7876.wiki.zoho.com/Introduction-to-R.html , а то у меня какое то чувство что я его зря переводил) Проблема в том, что это в среднем нужно рассчитать температуру за предыдущие семь дней, но это не всегда так. Есть периоды более 7 дней и менее 7 дней. Только для первого измерения всегда нужна средняя температура за предыдущие 7 дней. Меня то как раз рассчитать среднее за 7 дней очень устраивает:) Я не программист и в R только начала работать. Мои коллеги рекомендовали использовать for для того, чтобы рассчитать среднее... Почему вы пишите, что для for нужна крайне исключительная причина? Данные сейчас загружу |
|
29.08.2012 - 21:17
Сообщение
#6
|
|
Группа: Пользователи Сообщений: 43 Регистрация: 4.01.2012 Пользователь №: 23400 |
---
Сообщение отредактировал Olga_ - 8.11.2012 - 22:35 |
|
30.08.2012 - 10:44
Сообщение
#7
|
|
Группа: Пользователи Сообщений: 1091 Регистрация: 26.08.2010 Пользователь №: 22699 |
Проблема в том, что это в среднем нужно рассчитать температуру за предыдущие семь дней, но это не всегда так. Есть периоды более 7 дней и менее 7 дней. Только для первого измерения всегда нужна средняя температура за предыдущие 7 дней. Меня то как раз рассчитать среднее за 7 дней очень устраивает:) Я не программист и в R только начала работать. Мои коллеги рекомендовали использовать for для того, чтобы рассчитать среднее... Почему вы пишите, что для for нужна крайне исключительная причина? Данные сейчас загружу Данные это хорошо. Но: 1) Что значит "менее-более семи дней"?! Вы в первом сообщении написали а) что у Вас база данных температуры именно за _каждый_ день, б) что для каждого инциндента из таблицы случаев надо получить среднюю температуру за предыдущие 7 дней. 2) Объясните структуру своего эксперимента. Типа: "Для каждой строки из таблицы такой то получить то то и то то...". Ну или просто раскажите что надо получитью Иначе мне не понятно как надо группировать данные, и я не смогу Вам помочь. 3) Привет коллегам, они хорошие "программисты на фортране" В 95% случаев анализа данных в R for() не нужен (а следовательно , по заветам медицинской статистики не существует , именно для этого R и разрабатывали собственно --- что бы избежать низкоуровневого программирования. |
|
30.08.2012 - 12:02
Сообщение
#8
|
|
Группа: Пользователи Сообщений: 43 Регистрация: 4.01.2012 Пользователь №: 23400 |
Данные это хорошо. Но: 1) Что значит "менее-более семи дней"?! Вы в первом сообщении написали а) что у Вас база данных температуры именно за _каждый_ день, б) что для каждого инциндента из таблицы случаев надо получить среднюю температуру за предыдущие 7 дней. 2) Объясните структуру своего эксперимента. Типа: "Для каждой строки из таблицы такой то получить то то и то то...". Ну или просто раскажите что надо получитью Иначе мне не понятно как надо группировать данные, и я не смогу Вам помочь. 3) Привет коллегам, они хорошие "программисты на фортране" В 95% случаев анализа данных в R for() не нужен (а следовательно , по заветам медицинской статистики не существует , именно для этого R и разрабатывали собственно --- что бы избежать низкоуровневого программирования. Используя базу данных температуры за каждый день. (файл temp.txt ) нужно рассчитать среднюю температуру для каждого участника исследования (data.txt) за период с момента предыдущего наблюдения до следующего. Обычно разница между наблюдениями составляет 7 дней, но может быть и больше, например, для 1 пациента период между первым наблюдением и вторым stday 28apr1999 07may1999 Конечно, такие случаи встречаются редко и их можно просто игнорировать, используя предложенную вами методику рассчета (средняя температура за предыдущие 7 дней). Для первого измерения всегда нужна средняя температура за предыдущие 7 дней. Одна из целей исследования - оценить влияние факторов окружающей среды (темп., РМ10, О3 и т.д.) на частоту респираторных симптомов в течении 1 года. |
|
30.08.2012 - 15:37
Сообщение
#9
|
|
Группа: Пользователи Сообщений: 1091 Регистрация: 26.08.2010 Пользователь №: 22699 |
Используя базу данных температуры за каждый день. (файл temp.txt ) нужно рассчитать среднюю температуру для каждого участника исследования (data.txt) за период с момента предыдущего наблюдения до следующего. Обычно разница между наблюдениями составляет 7 дней, но может быть и больше, например, для 1 пациента период между первым наблюдением и вторым stday 28apr1999 07may1999 Конечно, такие случаи встречаются редко и их можно просто игнорировать, используя предложенную вами методику рассчета (средняя температура за предыдущие 7 дней). Для первого измерения всегда нужна средняя температура за предыдущие 7 дней. Одна из целей исследования - оценить влияние факторов окружающей среды (темп., РМ10, О3 и т.д.) на частоту респираторных симптомов в течении 1 года. Ольга, Вам не кажется что эти данные немного отличаются от первоначальных? Вы не ошиблись с файлом? Код day month year temp temp_min temp_max 1 1 1998 1.15 2 1 1998 6.64 3 1 1998 6.01 4 1 1998 6.16 5 1 1998 7.55 6 1 1998 6.36 7 1 1998 8.13 8 1 1998 7.66 9 1 1998 2.42 |
|
30.08.2012 - 16:10
Сообщение
#10
|
|
Группа: Пользователи Сообщений: 43 Регистрация: 4.01.2012 Пользователь №: 23400 |
Ольга, Вам не кажется что эти данные немного отличаются от первоначальных? Вы не ошиблись с файлом? Код day month year temp temp_min temp_max 1 1 1998 1.15 2 1 1998 6.64 3 1 1998 6.01 4 1 1998 6.16 5 1 1998 7.55 6 1 1998 6.36 7 1 1998 8.13 8 1 1998 7.66 9 1 1998 2.42 Нет, сначала , как пример, я привела модифицированные данные, где я обьединила в одну коллонку время исследования. |
|
30.08.2012 - 19:50
Сообщение
#11
|
|
Группа: Пользователи Сообщений: 1091 Регистрация: 26.08.2010 Пользователь №: 22699 |
Нет, сначала , как пример, я привела модифицированные данные, где я обьединила в одну коллонку время исследования. ну так почему бы и не выкладывать именно рабочий файл, раз уж Вы время привели к вменяемому формату? ну тогда уж продолжу прошлый "алгоритм": 1) получаем массив средних за 3-14 последних дней на текущую дату (колонку даты возьмете свою Код tmp<-sapply(1:10, function(i) avgtemp=runmean(temp.time$temp, i+3, alg= "exact", endrule="mean", align = "right") ) 2) разницу между текущим обследованием и предыдущим получаем Код > data<-read.table("data.txt", header=TRUE) > data$stday<-as.Date(data$stday, "%d%b%Y") > head(data) FID stday woche 1 1 1999-04-28 1 2 1 1999-05-07 2 3 1 1999-05-14 3 4 1 1999-05-21 4 5 1 1999-05-27 5 6 1 1999-06-04 6 > head(diff(data$stday)) Time differences in days [1] 9 7 7 6 8 7 Эти данные имеют на один день меньше, и их надо дополнить NA вначале если присоединять в data. 3) Берем из diff(data$stday)[i+1] где i номер строки из data имеющий смысл и получаем номер столбца скользящей средней температуры: tmp[,diff(data$stday)[i+1]] 4) Теперь вставляем выборку даты data$stday[i]==temp$day, получаем итоговую конструкцию tmp[ data$stday[i]==temp$day , diff(data$stday)[i+1] ] которая по i тому номеру строки из первой таблицы возвращает среднюю температуру за предыдущий период 5) осмысленные строки -- по красивому надо написать агрегирующую функцию в plyr но можно и проще Код sapply(2:nrow(data), function(i) if(data[i,1]==data[i-1,1]) tmp[ data$stday[i]==temp$day , diff(data$stday)[i+1] ] else NA )
Сообщение отредактировал p2004r - 30.08.2012 - 19:53 |
|
30.08.2012 - 21:16
Сообщение
#12
|
|
Группа: Пользователи Сообщений: 43 Регистрация: 4.01.2012 Пользователь №: 23400 |
ну так почему бы и не выкладывать именно рабочий файл, раз уж Вы время привели к вменяемому формату? ну тогда уж продолжу прошлый "алгоритм": 1) получаем массив средних за 3-14 последних дней на текущую дату (колонку даты возьмете свою Код tmp<-sapply(1:10, function(i) avgtemp=runmean(temp.time$temp, i+3, alg= "exact", endrule="mean", align = "right") ) 2) разницу между текущим обследованием и предыдущим получаем Код > data<-read.table("data.txt", header=TRUE) > data$stday<-as.Date(data$stday, "%d%b%Y") > head(data) FID stday woche 1 1 1999-04-28 1 2 1 1999-05-07 2 3 1 1999-05-14 3 4 1 1999-05-21 4 5 1 1999-05-27 5 6 1 1999-06-04 6 > head(diff(data$stday)) Time differences in days [1] 9 7 7 6 8 7 Эти данные имеют на один день меньше, и их надо дополнить NA вначале если присоединять в data. 3) Берем из diff(data$stday)[i+1] где i номер строки из data имеющий смысл и получаем номер столбца скользящей средней температуры: tmp[,diff(data$stday)[i+1]] 4) Теперь вставляем выборку даты data$stday[i]==temp$day, получаем итоговую конструкцию tmp[ data$stday[i]==temp$day , diff(data$stday)[i+1] ] которая по i тому номеру строки из первой таблицы возвращает среднюю температуру за предыдущий период 5) осмысленные строки -- по красивому надо написать агрегирующую функцию в plyr но можно и проще Код sapply(2:nrow(data), function(i) if(data[i,1]==data[i-1,1]) tmp[ data$stday[i]==temp$day , diff(data$stday)[i+1] ] else NA ) Огромное спасибо! То что я хотела, но никак не могла выразить из-за недостатка знаний в R. Вы, конечно, правы, нужно сначала читать, а потом еще перечитывать. Мне кажется, что вы меня уже не в первый раз выручаете. Как то регестрировалась я на одном форуме биологов н- ное время назад... Вчера файлы никак загружались, выложила, те файлы, которые удалось загрузить. Сообщение отредактировал Olga_ - 30.08.2012 - 21:17 |
|
11.09.2012 - 11:50
Сообщение
#13
|
|
Группа: Пользователи Сообщений: 43 Регистрация: 4.01.2012 Пользователь №: 23400 |
Что то у меня не получается..
Мой код airpollution <- read.table("...Payerne.csv",blank.lines.skip = TRUE,fill=T,header=T,sep=";") summary(airpollution) day <- paste(airpollution[,3],"-",airpollution[,2],"-",airpollution[,1],sep="") day <- as.Date(day, format="%Y-%m-%d") payerne <- cbind(airpollution,day) tmp<-sapply(1:10, function(i) avgtemp=runmean(payerne$temp, i+3, alg= "exact", endrule="mean", align = "right") ) data<- read.dta("...symptoms.dta") data$stday<-as.Date(data$stday, "%d%b%Y") head(diff(data$stday)) sapply(2:nrow(data), function(i) if(data[i,1]==data[i-1,1]) tmp[ data$stday[i]==payerne$day , diff(data$stday)[i+1] ] else NA ) Error in FUN(2:20426[[177L]], ...) : subscript out of bounds Сообщение отредактировал Olga_ - 11.09.2012 - 11:51 |
|
12.09.2012 - 12:30
Сообщение
#14
|
|
Группа: Пользователи Сообщений: 1091 Регистрация: 26.08.2010 Пользователь №: 22699 |
ну а файлы данных запаковать в rar и выложить?
|
|
6.11.2012 - 18:08
Сообщение
#15
|
|
Группа: Пользователи Сообщений: 43 Регистрация: 4.01.2012 Пользователь №: 23400 |
Решили оставить вариант с уровнем полютантов за предыдущие 7 дней.
Хотела бы использовать также для анализа время, в моем исследовании это 1 неделя. Но как было отмечено выше, время Исследований имеет неодинаковый промежуток. 1. Сначала, хотела бы создать одинаковые промежутки времени- недели начиная, к примеру с 01.01.1999. Какой это пакет в R? Где то кто то уже рекомендовал этот пакет, но не могу вспомнить 2. Затем интерполировать реальные даты в искусственные промежутки времени. С помощью какого пакета/пакетов R можно это выполнить? Спасибо заранее. Сообщение отредактировал Olga_ - 6.11.2012 - 18:09 |
|