Экстраполяция в R без использования Hmisc

Ответ на этот предыдущий вопрос: экстраполировать в R для данных временного ряда у меня не работает из-за R версий.

У меня есть фрейм данных NEI_othertier1_long, который выглядит примерно так:

state    pollutant    Sector       Fuel      description    year     value
AK       Ammonia      Refining     Diesel     industrial    2008      1.18
AK       Ammonia      Refining     Diesel     industrial    2009      NA
AK       Ammonia      Refining     Diesel     industrial    2010      NA
AK       Ammonia      Refining     Diesel     industrial    2011      5.76
AK       Ammonia      Refining     Diesel     industrial    2012      NA
AK       Ammonia      Refining     Diesel     industrial    2013      NA
AK       Ammonia      Refining     Diesel     industrial    2014      5.83
AK       Ammonia      Refining     Diesel     industrial    2015      NA
AK       Ammonia      Refining     Diesel     industrial    2016      NA
AK       Ammonia      Refining     Diesel     industrial    2017      8.96
AK       Ammonia      Refining     Diesel     industrial    2018      NA
AK       Ammonia      Refining     Diesel     industrial    2019      NA

У меня есть значения для 2008, 2011, 2014 и 2017 годов. Я смог успешно линейно интерполировать 2009-2016 годы, используя этот код:

    NEI_othertier1_long %>%
  dplyr::mutate( value = na.approx(value, na.rm = FALSE, rule = 2) ) -> NEI_othertier1_interpolated

Но интерполяция переносит значение 2017 года на 2018 и 2019 годы. Я хочу линейно экстраполировать значения 2018 и 2019 годов по сравнению с предыдущими годами.

У меня R версии 3.5.2 (и я не могу обновить), поэтому не могу установить latticeExtra, от которого Hmisc зависит использование функции approxExtrap.

Любая помощь приветствуется!

dput (head (NEI_othertier1_long)) структура (list (state = c (AK, AK, AK, AK, AK, AK, AK, AK, AK, AK, AK, AK), загрязнитель = c (Аммиак, Аммиак, Аммиак, Аммиак, Аммиак, Аммиак, Аммиак, Аммиак, Аммиак, Аммиак, Аммиак, Аммиак), CEDS_Sector = c (1A1b_Pet-refining, 1A1b_Pet-refining, 1A1b_Pet-refining, 1A1b_Pet-refining, 1A1b_Pet-refining, 1A1b_Pet-refining, refining 1A1b_Pet-refining, , 1A1b_Pet-refining, 1A1b_Pet-refining, 1A1b_Pet-refining, 1A1b_Pet-refining, 1A1b_Pet-refining), CEDS_Fuel = c (дизельное масло, дизельное масло, дизельное масло, дизельное масло, дизельное масло, дизельное масло, дизельное масло, дизельное масло, дизельное масло), дизельное масло , tier1_description = c (FUEL COMB. INDUSTRIAL, FUEL COMB. INDUSTRIAL, FUEL COMB. INDUSTRIAL, FUEL COMB. INDUSTRIAL, FUEL COMB. INDUSTRIAL, FUEL COMB. INDUSTRIAL, FUEL COMB. INDUSTRIAL, FUEL COMB. INDUSTRIAL. FUEL COMB. INDUSTRIAL, FUEL COMB. INDUSTRIAL, FUEL COMB. INDUSTRIAL), единица = c (TON, TON, TON, TON, TON, TON, TON, TON, TON, TON, TON, TON ), год = 2008: 2019, выбросы = c (1.18, NA, NA, 5.76, NA, NA, 5.83, NA, NA, 8.96, NA, NA)), row.names = c (NA, -12L), class = c (grouped_df, tbl_df, tbl, data.frame), groups = structure (list (state = AK, pollutant = Ammonia, CEDS_Sector = 1A1b_Pet-refining, CEDS_Fuel = diesel_oil, tier1_description = FUEL COMB. ПРОМЫШЛЕННЫЙ, unit = TON, .rows = list (1:12)), row.names = c (NA, -1L), class = c (tbl_df, tbl, data.frame), .drop = TRUE))


person Maridee Weber    schedule 14.09.2020    source источник


Ответы (1)


approxExtrap - это просто оболочка вокруг approx, поэтому вы можете просто скопировать определение функции и использовать его.

NEI_othertier1_long %>% dplyr::mutate(x = approxExtrap(year, value, year, na.rm = TRUE)$y)

Вот approxExtrap, если вы не можете его найти:

approxExtrap <- function (x, y, xout, method = "linear", n = 50, rule = 2, f = 0, 
  ties = "ordered", na.rm = FALSE) 
{
  if (is.list(x)) {
    y <- x[[2]]
    x <- x[[1]]
  }
  if (na.rm) {
    d <- !is.na(x + y)
    x <- x[d]
    y <- y[d]
  }
  d <- !duplicated(x)
  x <- x[d]
  y <- y[d]
  d <- order(x)
  x <- x[d]
  y <- y[d]
  w <- approx(x, y, xout = xout, method = method, n = n, rule = 2, 
    f = f, ties = ties)$y
  r <- range(x)
  d <- xout < r[1]
  if (any(is.na(d))) 
    stop("NAs not allowed in xout")
  if (any(d)) 
    w[d] <- (y[2] - y[1])/(x[2] - x[1]) * (xout[d] - x[1]) + 
    y[1]
  d <- xout > r[2]
  n <- length(y)
  if (any(d)) 
    w[d] <- (y[n] - y[n - 1])/(x[n] - x[n - 1]) * (xout[d] - 
        x[n - 1]) + y[n - 1]
  list(x = xout, y = w)
}
person pseudospin    schedule 14.09.2020
comment
Большое спасибо! Никогда бы не подумал об этом. - person Maridee Weber; 14.09.2020