Как найти среднее количество дней подряд, используя даты в R?
Я работаю с климатическими данными, к которым я добавил даты и разделил их на части, включив в них только данные, температура которых превышает 90-й процентиль в течение 3 или более дней.
set.seed(12368)
A <- data.frame("Y" = rep(1990:2009, each = 360), "M" = rep(1:12, each = 30, times = 20), "D" = rep(1:30, 240), "Temp" = c(sample(0:35, size = 7200, replace = TRUE)))
B <- data.frame("Y" = rep(1990:2009, each = 360), "M" = rep(1:12, each = 30, times = 20), "D" = rep(1:30, 240), "Temp" = c(sample(-10:22, size = 7200, replace = TRUE)))
C <- data.frame("Y" = rep(1990:2009, each = 360), "M" = rep(1:12, each = 30, times = 20), "D" = rep(1:30, 240), "Temp" = c(sample(3:42, size = 7200, replace = TRUE)))
climate <- list("Alist" = A, "Blist" = B, "Clist" = C)
#climate
library(dplyr)
library(magrittr)
datedclimate <- lapply(lapply(climate,
function(x)
x %<>%
mutate("Date" = as.Date(with(x, paste(Y,M,D, sep = "-")),"%Y-%m-%d"))
),
function(y)
y %>% relocate("Date")
)
datedclimate
tm <- lapply(datedclimate, \(x) {
y <- as.data.frame(subset(x,
Temp > quantile(Temp, probs = 0.90, na.rm = TRUE))
)
y[unique(
sort(
unlist(
lapply( # this iterates through all of my data sets
which(
c(diff(y[,'Date'])==1, FALSE) & c(diff(y[,'Date'], diff=2)==0, FALSE, FALSE)),
\(x) x + 0:2)
))),]
}
)
tm
Следующий шаг, который я собираюсь сделать, — это найти среднее количество дней для каждого последовательного количества дней. Например, подмножество списка C имеет 5 участков по 3 дня подряд и 1 участок по 4 дня подряд. Следовательно, среднее количество дней подряд составит 3,17 из (3+3+3+3+3+4)/6.
Я пытался изменить код из этого вопроса, но он не возвращает правильные числа. Я также предпочел бы изменить/поработать с подмножеством данных tm сверху, если это возможно.
ConsecMean <- function(x) {
x <- ifelse(x > quantile(x, probs = 0.90, na.rm = TRUE), 0, 1)
cs <- cumsum(x)
cs <- cs[x == 0]
mean <- mean(table(cs))
return(mean)
}
tri <- lapply(lapply(datedclimate, "[[", 5), ConsecMean)
tri
> tri
$Alist
[1] 1.105735
$Blist
[1] 1.104746
$Clist
[1] 1.099693
Итак, я ищу функцию, которую я могу применить к списку кадров данных (климат), которая будет возвращать это среднее число для каждого кадра данных в списке. Я думаю, что мне следует использовать Lapply, но я не знаю, что делать дальше. Я также пытался использовать rle, но он несовместим с форматом даты, и я также думаю, что rle — неподходящая функция для этого, лол. Какой код/функции вы бы порекомендовали использовать для решения этой проблемы?
Ответ всегда: используйте rle. Остальное — просто перестановка данных :-)
Ответы 3
Я думаю, что вы на правильном пути с rle(), но вам нужно использовать другую его версию. Если вы готовы рассмотреть подход на основе tibble, в котором используется групповая обработка, то я считаю, что мой код ниже сделает то, что вы хотите. Однако я не уверен в этом, поскольку предоставленный вами код генерации данных не соответствует вашему описанию того, что должно быть в Alist, когда я запускаю его на своем компьютере (R 4.4.1, Win 10).
Вы заявили:
подмножество списка A состоит из 4 участков по 3 дня подряд и 1 участка по 6 дней подряд. Следовательно, среднее количество дней подряд будет 3,6 из (3+3+3+3+6)/5.
Однако, когда я запускаю ваш код, я получаю шесть запусков в течение трех дней подряд; см. следующие (аннотированные) данные из tm$Alist:
> tm$Alist
Date Y M D Temp
2012 1995-08-02 1995 8 2 33 # start 1
2013 1995-08-03 1995 8 3 34
2014 1995-08-04 1995 8 4 33 # end 1
3641 2000-02-11 2000 2 11 34 # start 2
3642 2000-02-12 2000 2 12 33
3643 2000-02-13 2000 2 13 33 # end 2
3650 2000-02-20 2000 2 20 34 # start 3
3651 2000-02-21 2000 2 21 33
3652 2000-02-22 2000 2 22 35 # end 3
4066 2001-04-16 2001 4 16 34 # start 4
4067 2001-04-17 2001 4 17 33
4068 2001-04-18 2001 4 18 33 # end 4
4582 2002-09-22 2002 9 22 35 # start 5
4583 2002-09-23 2002 9 23 35
4584 2002-09-24 2002 9 24 33 # end 5
6024 2006-09-24 2006 9 24 34 # start 6
6025 2006-09-25 2006 9 25 35
6026 2006-09-26 2006 9 26 35 # end 6
Учитывая результат, который я вижу, правильное среднее количество последовательных дней будет ровно 3 для Alist. Если вы не хотите сгруппировать данные по году и месяцу, в этом случае прогоны № 2 и 3, приведенные выше, свернутся в один прогон продолжительностью шесть дней?
В любом случае, вот код, который я придумал. Обратите внимание, что он использует функцию rleidv() из пакета data.table, но я совершенно уверен, что data.table устанавливается при установке пакета tidyverse. Другими словами, я предполагаю, что он у вас уже установлен.
set.seed(12368)
A <- data.frame("Y" = rep(1990:2009, each = 360), "M" = rep(1:12, each = 30, times = 20), "D" = rep(1:30, 240), "Temp" = c(sample(0:35, size = 7200, replace = TRUE)))
B <- data.frame("Y" = rep(1990:2009, each = 360), "M" = rep(1:12, each = 30, times = 20), "D" = rep(1:30, 240), "Temp" = c(sample(-10:22, size = 7200, replace = TRUE)))
C <- data.frame("Y" = rep(1990:2009, each = 360), "M" = rep(1:12, each = 30, times = 20), "D" = rep(1:30, 240), "Temp" = c(sample(3:42, size = 7200, replace = TRUE)))
library(data.table)
library(dplyr)
bind_rows(A, B, C, .id = "id") |>
as_tibble() |>
mutate("Date" = as.Date(paste(Y, M, D, sep = "-"), "%Y-%m-%d")) |>
filter(is.na(Date) == FALSE) |>
relocate(Date) |>
arrange(id, Date) |>
mutate(
keepflag = Temp > quantile(Temp, probs = 0.90, na.rm = TRUE),
runs = data.table::rleidv(keepflag),
.by = id
) |>
filter(keepflag == TRUE) |>
count(id, runs) |>
filter(n >= 3) |>
summarize(ConsecMean = mean(n), .by = id)
В приведенном выше коде:
bind_rows()«складывает» data.frames один поверх другого и создает дополнительную переменную с именемid, чтобы отслеживать, из какого data.frame взята каждая строка в выходных данных.as_tibble()преобразует data.frame в тиббл для более удобной печати- строит переменную
Dateиз компонентовY,MиD filter()сбрасывает записи, в которых построена датаNA, а именно февраль каждого года.relocate()перемещает переменнуюDateв первый столбец тибблаarrange()сначала сортирует тиббл по исходной таблице (id), затем поDateпри подготовке к использованию группового кодирования длин серий.- второй
mutate()обрабатывает данные по значению столбцаid(.by = id) для создания (1) логического флага, который идентифицирует записи, находящиеся выше 90-го процентиля дляTemp, и (2) идентификатора RLE (с использованием функцииdata.table::rleidv()), который присваивает целое число ( начиная с 1) для каждой записи, увеличиваясь на 1 при обнаружении другого значения, эффективно группируя последовательные даты filter(keepflag == TRUE)удаляет все записи, гдеTempне превышал 90-го процентиляcount()подсчитывает количество дат в каждом наборе идентификаторов RLE из каждого исходного data.frames и сохраняет это значение (по умолчанию) в новой переменной с именемnfilter(n >= 3)отбрасывает пробежки, длившиеся не менее трех дней подряд, иsummarize()создает переменнуюConsecMeanпо значению столбцаid
Результат, который я получаю для ваших данных, используя этот код:
# A tibble: 3 × 2
id ConsecMean
<chr> <dbl>
1 1 3
2 2 3
3 3 3.17
Привет, да, ты прав насчет списка А, я случайно его неправильно прочитал, моя вина.
Возможно, вы захотите использовать некоторые части этого подхода. Обратите внимание, что нет необходимости во внешних пакетах.
lapply(climate, \(l) {
x = transform(l, Date = as.Date(paste(Y, M, D, sep = "-"), "%Y-%m-%d")) |>
subset(select = c(Date, Temp), # cosmetics
subset = Temp > quantile(Temp, probs = .9, na.rm = TRUE)) |>
sort_by(~Date)
i = cumsum(c(0L, diff(x$Date) != 1L))
data.frame(start = x$Date[c(1L, diff(i)) == 1L], n = rle(i)$lengths) |>
subset(n > 2L) |> with(mean(n))
})
предоставление
$Alist
[1] 3
$Blist
[1] 3
$Clist
[1] 3.166667
Я обернул части обработки ненужных данных и некоторые косметические детали внутри одного lapply-вызова. Возможно, вы захотите использовать их в другом месте вашего кода, возможно, отдельно.
Если убрать хитрый |> with(mean(n)), возвращаются удобные фреймы данных, содержащие даты начала:
$Alist
start n
142 1995-08-02 3
271 2000-02-11 3
272 2000-02-20 3
306 2001-04-16 3
354 2002-09-22 3
466 2006-09-24 3
$Blist
start n
50 1992-01-07 3
204 1996-11-20 3
278 1999-06-17 3
413 2003-11-27 3
484 2005-11-29 3
$Clist
start n
25 1990-07-19 3
33 1990-11-08 3
266 1997-08-17 3
320 1999-02-13 3
325 1999-04-19 3
555 2007-03-19 4
Использование tm в примечании в конце seqid присваивает уникальный идентификатор, который одинаков для всех дат в любой конкретной последовательной последовательности дат, table подсчитывает элементы в каждой такой группе, давая вектор счетчиков, и mean берет их среднее значение.
library(collapse)
sapply(tm, function(x) mean(table(seqid(x$Date))))
## Alist Blist Clist
## 3.000000 3.000000 3.166667
Примечание
tm <- list(Alist = structure(list(Date = structure(c(9344, 9345, 9346,
10998, 10999, 11000, 11007, 11008, 11009, 11428, 11429, 11430,
11952, 11953, 11954, 13415, 13416, 13417), class = "Date"), Y = c(1995L,
1995L, 1995L, 2000L, 2000L, 2000L, 2000L, 2000L, 2000L, 2001L,
2001L, 2001L, 2002L, 2002L, 2002L, 2006L, 2006L, 2006L), M = c(8L,
8L, 8L, 2L, 2L, 2L, 2L, 2L, 2L, 4L, 4L, 4L, 9L, 9L, 9L, 9L, 9L,
9L), D = c(2L, 3L, 4L, 11L, 12L, 13L, 20L, 21L, 22L, 16L, 17L,
18L, 22L, 23L, 24L, 24L, 25L, 26L), Temp = c(33L, 34L, 33L, 34L,
33L, 33L, 34L, 33L, 35L, 34L, 33L, 33L, 35L, 35L, 33L, 34L, 35L,
35L)), row.names = c(2012L, 2013L, 2014L, 3641L, 3642L, 3643L,
3650L, 3651L, 3652L, 4066L, 4067L, 4068L, 4582L, 4583L, 4584L,
6024L, 6025L, 6026L), class = "data.frame"), Blist = structure(list(
Date = structure(c(8041, 8042, 8043, 9820, 9821, 9822, 10759,
10760, 10761, 12383, 12384, 12385, 13116, 13117, 13118), class = "Date"),
Y = c(1992L, 1992L, 1992L, 1996L, 1996L, 1996L, 1999L, 1999L,
1999L, 2003L, 2003L, 2003L, 2005L, 2005L, 2005L), M = c(1L,
1L, 1L, 11L, 11L, 11L, 6L, 6L, 6L, 11L, 11L, 11L, 11L, 11L,
12L), D = c(7L, 8L, 9L, 20L, 21L, 22L, 17L, 18L, 19L, 27L,
28L, 29L, 29L, 30L, 1L), Temp = c(20L, 20L, 22L, 22L, 20L,
20L, 20L, 20L, 21L, 20L, 20L, 20L, 20L, 21L, 21L)), row.names = c(727L,
728L, 729L, 2480L, 2481L, 2482L, 3407L, 3408L, 3409L, 5007L,
5008L, 5009L, 5729L, 5730L, 5731L), class = "data.frame"), Clist = structure(list(
Date = structure(c(7504, 7505, 7506, 7616, 7617, 7618, 10090,
10091, 10092, 10635, 10636, 10637, 10700, 10701, 10702, 13591,
13592, 13593, 13594), class = "Date"), Y = c(1990L, 1990L,
1990L, 1990L, 1990L, 1990L, 1997L, 1997L, 1997L, 1999L, 1999L,
1999L, 1999L, 1999L, 1999L, 2007L, 2007L, 2007L, 2007L),
M = c(7L, 7L, 7L, 11L, 11L, 11L, 8L, 8L, 8L, 2L, 2L, 2L,
4L, 4L, 4L, 3L, 3L, 3L, 3L), D = c(19L, 20L, 21L, 8L, 9L,
10L, 17L, 18L, 19L, 13L, 14L, 15L, 19L, 20L, 21L, 19L, 20L,
21L, 22L), Temp = c(39L, 39L, 40L, 39L, 40L, 41L, 39L, 39L,
42L, 41L, 41L, 42L, 40L, 40L, 42L, 42L, 42L, 40L, 41L)), row.names = c(199L,
200L, 201L, 308L, 309L, 310L, 2747L, 2748L, 2749L, 3283L, 3284L,
3285L, 3349L, 3350L, 3351L, 6199L, 6200L, 6201L, 6202L), class = "data.frame"))
Кроме того: последние две строки (
mean <- mean(table(cs)),return(mean)) содержат много синтаксического шума и беспорядка дляmean(table(cs))(нет необходимости ни во временной переменной , ни в return() ). Отсутствие ненужного шума в коде облегчает его чтение и поддержку, а также значительно снижает количество ошибок.