У мене є кадр даних, який містить близько 200стовпці, що представляють щомісячні вимірювання посухи за період з 1999 по 2015 рік. Значення в кожному стовпчику можуть бути як позитивними, так і негативними. Кожен рядок у кадрі даних являє собою фокусний рік, на який я зацікавлений у обчисленні показників посилань. Фокусний рік представлений в іншому стовпчику. Можливо, буде кілька рядків з тим самим фокусним роком, якщо вони представляють вимірювання з різних сайтів (стовпець FIPS Тут "версія іграшки (оновлена версія!) Фрейму даних:
structure(list(FIPS = c(19045, 48157, 20045, 20027), Year = c(2003,
2004, 2005, 2005), pdsi_2002.01.15 = c(1.46, 4.38, 0.38, -1.41
), pdsi_2002.02.15 = c(1.6, 3.63, -0.05, -1.66), pdsi_2002.03.15 = c(1.32,
3, -0.62, -1.93), pdsi_2002.04.15 = c(1.81, 2.68, 0.66, -1.88
), pdsi_2002.05.15 = c(2.03, 1.86, 1.26, -1.7), pdsi_2002.06.15 = c(2.51,
1.74, -0.5, -2.94), pdsi_2002.07.15 = c(2.79, 1.94, -1.47, -3.82
), pdsi_2002.08.15 = c(3.06, 2.64, -1.99, -4.09), pdsi_2002.09.15 = c(2.08,
3.02, -2.82, -4.87), pdsi_2002.10.15 = c(2.68, 4.73, -2.02, -3.01
), pdsi_2002.11.15 = c(2, 5.28, -2.55, -3.22), pdsi_2002.12.15 = c(1.55,
5.94, -3.23, -3.52), pdsi_2003.01.15 = c(0.96, 5.39, -3.58, -3.51
), pdsi_2003.02.15 = c(0.29, 5.24, -3.54, -3.29), pdsi_2003.03.15 = c(-0.15,
4.41, -3.77, -3.15), pdsi_2003.04.15 = c(-1.13, 3.39, -3.33,
-2.46), pdsi_2003.05.15 = c(-1.05, 1.91, -3.47, -2.63), pdsi_2003.06.15 = c(-1.5,
1.45, -2.94, -2.34), pdsi_2003.07.15 = c(-0.85, 1.69, -3.42,
-3.02), pdsi_2003.08.15 = c(-1.78, 1.48, -2.75, -3.13), pdsi_2003.09.15 = c(-1.55,
2.31, -2.66, -2.85), pdsi_2003.10.15 = c(-1.87, 2.5, -2.99, -3.16
), pdsi_2003.11.15 = c(-1.19, 2.72, -3.39, -2.73), pdsi_2003.12.15 = c(0.09,
2.67, -2.96, -2.63), pdsi_2004.01.15 = c(-0.2, 3.2, -2.83, -2.42
), pdsi_2004.02.15 = c(0.07, 3.73, -2.78, -2.21), pdsi_2004.03.15 = c(1.58,
3.04, -1.66, -0.77), pdsi_2004.04.15 = c(0.37, 3.19, -2, -1.25
), pdsi_2004.05.15 = c(1.7, 3.71, -1.35, -1.41), pdsi_2004.06.15 = c(1.53,
5.21, -0.84, -1.04), pdsi_2004.07.15 = c(1.14, 4.84, 2.08, 0.93
), pdsi_2004.08.15 = c(1.4, 4.41, 3.22, 0.24), pdsi_2004.09.15 = c(-0.43,
3.27, 2.39, -0.44), pdsi_2004.10.15 = c(0.77, 2.77, 2.49, -1.11
), pdsi_2004.11.15 = c(0.94, 4.95, 2.94, -1.03), pdsi_2004.12.15 = c(0.62,
4.41, 2.67, -1.43), pdsi_2005.01.15 = c(1.51, 3.93, 3.55, -1.05
), pdsi_2005.02.15 = c(1.45, 4.54, 3.83, 0.71), pdsi_2005.03.15 = c(0.58,
4.31, 3.01, 0.24), pdsi_2005.04.15 = c(-0.97, 3.36, 1.97, 0.94
), pdsi_2005.05.15 = c(-1.57, 3.12, 1.54, -0.33), pdsi_2005.06.15 = c(-2.65,
2.02, 2.33, 1.16), pdsi_2005.07.15 = c(-3.58, 2.07, 2.31, 1.08
), pdsi_2005.08.15 = c(-3.51, 1.56, 3.7, 1.72), pdsi_2005.09.15 = c(-3.96,
-0.71, 3.62, 0.74), pdsi_2005.10.15 = c(-4.77, -2.13, 3.79, 0.96
), pdsi_2005.11.15 = c(-5.08, -2.32, 3.4, 0.53), pdsi_2005.12.15 = c(-5.63,
-2.57, 3.27, -0.22)), .Names = c("FIPS", "Year", "pdsi_2002.01.15",
"pdsi_2002.02.15", "pdsi_2002.03.15", "pdsi_2002.04.15", "pdsi_2002.05.15",
"pdsi_2002.06.15", "pdsi_2002.07.15", "pdsi_2002.08.15", "pdsi_2002.09.15",
"pdsi_2002.10.15", "pdsi_2002.11.15", "pdsi_2002.12.15", "pdsi_2003.01.15",
"pdsi_2003.02.15", "pdsi_2003.03.15", "pdsi_2003.04.15", "pdsi_2003.05.15",
"pdsi_2003.06.15", "pdsi_2003.07.15", "pdsi_2003.08.15", "pdsi_2003.09.15",
"pdsi_2003.10.15", "pdsi_2003.11.15", "pdsi_2003.12.15", "pdsi_2004.01.15",
"pdsi_2004.02.15", "pdsi_2004.03.15", "pdsi_2004.04.15", "pdsi_2004.05.15",
"pdsi_2004.06.15", "pdsi_2004.07.15", "pdsi_2004.08.15", "pdsi_2004.09.15",
"pdsi_2004.10.15", "pdsi_2004.11.15", "pdsi_2004.12.15", "pdsi_2005.01.15",
"pdsi_2005.02.15", "pdsi_2005.03.15", "pdsi_2005.04.15", "pdsi_2005.05.15",
"pdsi_2005.06.15", "pdsi_2005.07.15", "pdsi_2005.08.15", "pdsi_2005.09.15",
"pdsi_2005.10.15", "pdsi_2005.11.15", "pdsi_2005.12.15"), row.names = c(13222L,
18125L, 19543L, 19534L), class = "data.frame")
Що я хотів би зробити, це обчислити довжину іпідсумок кожного прогону негативних значень у фокусному році (тому шукаємо прогони в тому ж рядку в стовпцях), потім обчислюємо середню довжину пробігу, середню суму пробігу та середню суму кожної пробігу, поділену на кожну довжину пробігу для кожного рядка . Якщо додати ще один рівень складності, якщо вимірювання січня у фокальному році є негативним, я "хотів би потім озирнутися до років перед фокальним роком, щоб пояснити ситуацію, коли кількість негативних чисел почалася в попередньому році. Пробіг міг би продовжити весь шлях до січня 1999 року.
Я "зміг обчислити метрику довжини пробігу за допомогою rle (), але не зміг зрозуміти, як отримати суму пробігу.
Відповіді:
1 для відповіді № 1Я думаю, що це може працювати для того, що ти шукаєшдля цього, це буде генерувати 3 необхідних значення для вказаного року, і якщо в січні буде негативне значення, воно продовжуватиме зменшуватися до досягнення позитивного значення в попередньому році.
library(tidyr)
library(dplyr)
select.order <- colnames(drought_data)[3:length(colnames(drought_data))]
drought_data <- drought_data %>%
# Gather data by date
gather(key = date, value = value, -Year, -FIPS) %>%
# Separate date into separate columns
separate(date, into = c("yr","month", "day"), sep = "\.") %>%
# Extract year
mutate(yr = substr(yr, 6, 9)) %>%
# Sort data by FIPS number, year, month
arrange(FIPS, yr, month) %>%
# Group data by FIPS number, focal year, and data year
group_by(FIPS, Year, yr) %>%
# Generate a run number for each run of negative numbers for the focal year
mutate(run.num = ifelse(Year == yr,
{run.num = rle(ifelse(value < 0, 1, 0))
rep(ifelse(run.num$values == 1, cumsum(run.num$values), 0), run.num$lengths)}, NA),
# Set run.num to -1 for positive values
run.num = ifelse(value >= 0, -1, run.num)) %>%
# Sort data by FIPS number, descending year, and descending month
arrange(FIPS, desc(yr), desc(month)) %>%
# Group data by FIPS number and focal year
group_by(FIPS, Year) %>%
# Fill out the run numbers for each run to cross data years
fill(run.num, .direction = "down") %>%
# Convert all -1 run numbers (Which indicate positive values) to zero
mutate(run.num = ifelse(run.num == -1, 0, run.num),
# Set run.num for negative values that did not qualify as a run for the specified year to 0
run.num = ifelse(is.na(run.num), 0, run.num)) %>%
ungroup %>%
# mutate(run.num = ifelse(is.na(run.num, 0, run.num))) %>%
# Group data by FIPS number, focal year, and run number
group_by(FIPS, Year, run.num) %>%
# Calculate the length, sum, and rate of each run
mutate(run.length = ifelse(run.num == 0, 0, n()),
run.sum = ifelse(run.num == 0, 0, sum(value)),
run.rate = ifelse(run.num == 0, 0, run.sum/run.length)) %>%
# Group by FIPS number and focal year
group_by(FIPS, Year) %>%
# Calculate the mean run length, and mean run sum for the focal year of each FIPS number
mutate(mean.run.length = sum(ifelse(run.num == 0, 0, 1)) / max(run.num),
mean.run.length = ifelse(is.nan(mean.run.length), 0, mean.run.length),
mean.run.sum = sum(ifelse(run.num == 0, 0, value) / max(run.num)),
mean.run.sum = ifelse(is.nan(mean.run.sum), 0, mean.run.sum)) %>%
# Combine date parts back to single column
unite(dt, yr:day, sep = ".") %>%
# Recreate the pdsi_ label format on the date column
mutate(dt = paste0("pdsi_", dt)) %>%
# Drop the run.sum column
select(-run.sum) %>%
# Spread the data back to a wide view to eliminate duplicate run.rate values
spread(dt, value) %>%
# Group data by FIPS number and focal year
group_by(FIPS, Year) %>%
# Calculate the mean of the sum of run rates over the number of runs
mutate(mean.run.sum.length = sum(run.rate) / max(run.num),
mean.run.sum.length = ifelse(is.nan(mean.run.sum.length), 0, mean.run.sum.length)) %>%
# Remove grouping
ungroup %>%
# Drop the run.num, run.length, and run.rate columns
select(-run.num, -run.length, -run.rate) %>%
# Gather the data into tall view to remove duplicates and NA values
gather_("dt", "value", select.order, na.rm = TRUE) %>%
# Spread data back to wide view
spread(dt, value)
# Change the column order
drought_data <- drought_data[,c("FIPS","Year","mean.run.length","mean.run.sum","mean.run.sum.length", select.order)]
Кінцевим результатом буде оригінальний кадр даних із включеними трьома додатковими обчисленими стовпцями. Нижче наведено висновок обчислених стовпців для наданого тестового набору даних.
> drought_data[,c("FIPS","Year","mean.run.length","mean.run.sum","mean.run.sum.length")]
# A tibble: 4 x 5
FIPS Year mean.run.length mean.run.sum mean.run.sum.length
<dbl> <dbl> <dbl> <dbl> <dbl>
1 19045 2003 9.000000 -11.07 -1.2300000
2 20027 2005 2.333333 -1.87 -0.5206667
3 20045 2005 0.000000 0.00 0.0000000
4 48157 2004 0.000000 0.00 0.0000000
0 для відповіді № 2
Ось а tidyverse підхід до проблеми, принаймні перша частина. Але я вважаю, що це потрапляє і у другу частину вашої проблеми.
На мою думку, корисно привести дані в іншому, охайному форматі, де кожен рядок - щомісячне спостереження в окремому фокусному році. (До речі, я не був впевнений, чому focal_ years колона була вся 2001. Згідно з вашим описом даних, чи не повинні бути окремі роки?)
library(tidyverse)
drought_data_tidy <- drought_data %>%
gather(key, value, -FIPS, -Year) %>%
arrange(FIPS)
Це дає нам таку структуру даних (зауважте, що я додав ідентифікатор рядка, оскільки фокусний рік завжди однаковий. Якщо це помилка з вашого боку, код дещо спрощує):
> head(drought_data_tidy)
FIPS Year key value
1 8019 2005 pdsi_2002.01.15 -1.73
2 8019 2005 pdsi_2002.02.15 -2.04
3 8019 2005 pdsi_2002.03.15 -2.44
4 8019 2005 pdsi_2002.04.15 -3.55
5 8019 2005 pdsi_2002.05.15 -3.84
6 8019 2005 pdsi_2002.06.15 -4.42
Єдина «складність» - створити унікальний run_id для кожного негативного циклу, який ми отримуємо таким чином:
drought_data_tidy <- drought_data_tidy %>%
group_by(FIPS) %>%
mutate(run_id = cumsum(c(TRUE, diff(value < 0) != 0)),
run_id = ifelse(value < 0, run_id, NA))
Тепер все, що залишається, є відповідним group_by і mean команди. Зауважте, що я вибираю відповідні частини значень, за якими ми робимо обчислення [!is.na(run_id)]; ми також могли б створити для цього ще одну змінну групування, яка може бути більш елегантною.
З огляду на ваші оновлені дані та проблеми з unique() Я вирішив провести обчислення summarize() команди та зберігати результати в окремих данихкадри, які я приєдную до кінця. Це також можливо зробити в одному гігантському ланцюжку dplyr, але я думаю, що цей підхід дещо простіший для читання і краще для виправлення помилок.
# run length
drought_data_run_length <- drought_data_tidy %>%
group_by(FIPS, run_id) %>%
summarize(run_length = n()) %>%
mutate(mean_run_length = mean(run_length[!is.na(run_id)]))
# mean run length for join
drought_data_mean_run_length <- drought_data_run_length %>%
group_by(FIPS) %>%
summarise(mean_run_length = unique(mean_run_length))
# run sum
drought_data_tidy <- drought_data_tidy %>%
group_by(FIPS, run_id) %>%
mutate(run_sum = sum(value))
# mean run sum
drought_data_mean_run_sum <- drought_data_tidy %>%
group_by(FIPS) %>%
summarise(mean_run_sum = mean(run_sum[!is.na(run_id)]))
# mean run sum by mean run length
drought_data_mrs_by_mrl <- left_join(drought_data_mean_run_sum,
drought_data_mean_run_length,
by = "FIPS") %>%
mutate(mrs_by_mrl = mean(mean_run_sum / mean_run_length))
# join run length, mean run length, mean run sum, mrs_by_mrl
drought_data_tidy <- left_join(drought_data_tidy,
drought_data_run_length %>% select(-mean_run_length),
by = c("FIPS", "run_id"))
drought_data_tidy <- left_join(drought_data_tidy,
drought_data_mean_run_length %>% select(FIPS, mean_run_length),
by = "FIPS")
drought_data_tidy <- left_join(drought_data_tidy,
drought_data_mean_run_sum %>% select(FIPS, mean_run_sum),
by = "FIPS")
drought_data_tidy <- left_join(drought_data_tidy,
drought_data_mrs_by_mrl %>% select(FIPS, mrs_by_mrl),
by = "FIPS")
Це дає наступний вихід:
> head(drought_data_tidy)
Source: local data frame [6 x 10]
Groups: FIPS, run_id [1]
FIPS Year key value run_id run_sum run_length mean_run_length mean_run_sum mrs_by_mrl
<dbl> <dbl> <chr> <dbl> <int> <dbl> <int> <dbl> <dbl> <dbl>
1 8019 2005 pdsi_2002.01.15 -1.73 1 -49.33 14 9 -33.54481 -3.542602
2 8019 2005 pdsi_2002.02.15 -2.04 1 -49.33 14 9 -33.54481 -3.542602
3 8019 2005 pdsi_2002.03.15 -2.44 1 -49.33 14 9 -33.54481 -3.542602
4 8019 2005 pdsi_2002.04.15 -3.55 1 -49.33 14 9 -33.54481 -3.542602
5 8019 2005 pdsi_2002.05.15 -3.84 1 -49.33 14 9 -33.54481 -3.542602
6 8019 2005 pdsi_2002.06.15 -4.42 1 -49.33 14 9 -33.54481 -3.542602
Тепер ви можете легко поширити це назад у вихідний формат даних. Однак я вважаю, що tidy Кадр даних щомісячних спостережень у довгому форматі корисніше отримати у другій частині вашого запитання з розумним використанням first(). Але для вирішення цього нам знадобиться трохи більше ваших даних (або даних про іграшки), які поширюються на різні focal_years.
Надія, що допомагає.