data.table линейно интерполируя значения NA без групп

18

Я хотел заполнить некоторые значения NA в data.table без групп. Пожалуйста, рассмотрите этот фрагмент data.table, представляющий время и расстояния:

library(data.table)
df <- data.frame(time = seq(7173, 7195, 1), dist = c(31091.33, NA, 31100.00, 31103.27, NA, NA, NA, NA, 31124.98, NA,31132.81, NA, NA, NA, NA, 31154.19, NA, 31161.47, NA, NA, NA, NA, 31182.97))
DT<- data.table(df)

Я хочу в DT data.table заполнить значения NA функцией, зависящей от значения не-NA до / после. Например, написание функции в j для замены каждой инструкции

DT[2, dist := (31091.33 + (31100-31091.33) / 2)]

тогда

DT[5:8, dist := (31103.27 + "something" * (31124.98 - 31103.27) / 5)]

так далее...

r data.table ArnaudR
источник

7

Код объясняется в строке. Вы можете удалить временные столбцы, используя df[,dist_before := NULL], например.

library(data.table)
df=data.table(time=seq(7173,7195,1),dist=c(31091.33,NA,31100.00,31103.27,NA,NA,NA,
NA,31124.98,NA,31132.81,NA,NA,NA,NA,31154.19,NA,31161.47,NA,NA,NA,NA,31182.97))
df
#>     time     dist
#>  1: 7173 31091.33
#>  2: 7174       NA
#>  3: 7175 31100.00
#>  4: 7176 31103.27
#>  5: 7177       NA
#>  6: 7178       NA
#>  7: 7179       NA
#>  8: 7180       NA
#>  9: 7181 31124.98
#> 10: 7182       NA
#> 11: 7183 31132.81
#> 12: 7184       NA
#> 13: 7185       NA
#> 14: 7186       NA
#> 15: 7187       NA
#> 16: 7188 31154.19
#> 17: 7189       NA
#> 18: 7190 31161.47
#> 19: 7191       NA
#> 20: 7192       NA
#> 21: 7193       NA
#> 22: 7194       NA
#> 23: 7195 31182.97
#>     time     dist
# Carry forward the last non-missing observation
df[,dist_before := nafill(dist, "locf")]
# Bring back the next non-missing dist
df[,dist_after := nafill(dist, "nocb")]
# rleid will create groups based on run-lengths of values within the data.
# This means 4 NA's in a row will be grouped together, for example.
# We then count the missings and add 1, because we want the 
# last NA before the next non-missing to be less than the non-missing value.
df[, rle := rleid(dist)][,missings := max(.N +  1 , 2), by = rle][]
#>     time     dist dist_before dist_after rle missings
#>  1: 7173 31091.33    31091.33   31091.33   1        2
#>  2: 7174       NA    31091.33   31100.00   2        2
#>  3: 7175 31100.00    31100.00   31100.00   3        2
#>  4: 7176 31103.27    31103.27   31103.27   4        2
#>  5: 7177       NA    31103.27   31124.98   5        5
#>  6: 7178       NA    31103.27   31124.98   5        5
#>  7: 7179       NA    31103.27   31124.98   5        5
#>  8: 7180       NA    31103.27   31124.98   5        5
#>  9: 7181 31124.98    31124.98   31124.98   6        2
#> 10: 7182       NA    31124.98   31132.81   7        2
#> 11: 7183 31132.81    31132.81   31132.81   8        2
#> 12: 7184       NA    31132.81   31154.19   9        5
#> 13: 7185       NA    31132.81   31154.19   9        5
#> 14: 7186       NA    31132.81   31154.19   9        5
#> 15: 7187       NA    31132.81   31154.19   9        5
#> 16: 7188 31154.19    31154.19   31154.19  10        2
#> 17: 7189       NA    31154.19   31161.47  11        2
#> 18: 7190 31161.47    31161.47   31161.47  12        2
#> 19: 7191       NA    31161.47   31182.97  13        5
#> 20: 7192       NA    31161.47   31182.97  13        5
#> 21: 7193       NA    31161.47   31182.97  13        5
#> 22: 7194       NA    31161.47   31182.97  13        5
#> 23: 7195 31182.97    31182.97   31182.97  14        2
#>     time     dist dist_before dist_after rle missings
# .SD[,.I] will get us the row number relative to the group it is in. 
# For example, row 5 dist is calculated as
# dist_before + 1 * (dist_after - dist_before)/5
df[is.na(dist), dist := dist_before + .SD[,.I] *
                     (dist_after - dist_before)/(missings), by = rle]
df[]
#>     time     dist dist_before dist_after rle missings
#>  1: 7173 31091.33    31091.33   31091.33   1        2
#>  2: 7174 31095.67    31091.33   31100.00   2        2
#>  3: 7175 31100.00    31100.00   31100.00   3        2
#>  4: 7176 31103.27    31103.27   31103.27   4        2
#>  5: 7177 31107.61    31103.27   31124.98   5        5
#>  6: 7178 31111.95    31103.27   31124.98   5        5
#>  7: 7179 31116.30    31103.27   31124.98   5        5
#>  8: 7180 31120.64    31103.27   31124.98   5        5
#>  9: 7181 31124.98    31124.98   31124.98   6        2
#> 10: 7182 31128.90    31124.98   31132.81   7        2
#> 11: 7183 31132.81    31132.81   31132.81   8        2
#> 12: 7184 31137.09    31132.81   31154.19   9        5
#> 13: 7185 31141.36    31132.81   31154.19   9        5
#> 14: 7186 31145.64    31132.81   31154.19   9        5
#> 15: 7187 31149.91    31132.81   31154.19   9        5
#> 16: 7188 31154.19    31154.19   31154.19  10        2
#> 17: 7189 31157.83    31154.19   31161.47  11        2
#> 18: 7190 31161.47    31161.47   31161.47  12        2
#> 19: 7191 31165.77    31161.47   31182.97  13        5
#> 20: 7192 31170.07    31161.47   31182.97  13        5
#> 21: 7193 31174.37    31161.47   31182.97  13        5
#> 22: 7194 31178.67    31161.47   31182.97  13        5
#> 23: 7195 31182.97    31182.97   31182.97  14        2
#>     time     dist dist_before dist_after rle missings

smingerson
источник

8

Вы можете использовать approxфункцию для линейной интерполяции.

Для каждой группы NAs получите это подмножество DTплюс строки до и после. Затем примените approxк этому подмножеству distвектора, nаргумент которого approxравен числу строк в подмножестве .N.

DT[, g := rleid(dist)]

DT[is.na(dist), dist := {
      i <- .I[c(1, .N)] + c(-1, 1)
      DT[i[1]:i[2], approx(dist, n = .N)$y[-c(1, .N)]]
  }, by = g]

Или без approx

DT[, g := rleid(dist)]

DT[is.na(dist), dist := {
      i <- .I[c(1, .N)] + c(-1, 1)
      DT[i[1]:i[2], dist[1] + 1:(.N - 2)*(dist[.N] - dist[1])/(.N - 1)]
  }, by = g]

редактировать: так как этот ответ был принят, я чувствую, что должен отметить, что другие ответы быстрее, и вторая часть ответа @ dww - это, по сути, мой первый блок кода, но с удаленной ненужной частью группировки (так что это проще и быстрее).

IceCreamToucan
источник

На самом деле, я задаю этот вопрос, а потом пытаюсь сделать нелинейное приближение, чтобы ваше решение было более приспособленным к моим потребностям. вот почему я принял ваше решение

ArnaudR

6

2 других варианта:

1) подвижное соединение:

DT[is.na(dist), dist := {
        x0y0 <- DT[!is.na(dist)][.SD, on=.(time), roll=Inf, .(time=x.time, dist=x.dist)]
        x1y1 <- DT[!is.na(dist)][.SD, on=.(time), roll=-Inf, .(time=x.time, dist=x.dist)]
        (x1y1$dist - x0y0$dist) / (x1y1$time - x0y0$time) * (time - x0y0$time) + x0y0$dist
    }]
DT

2) другой близкий вариант ответа smingerson с использованием nafill

DT[, dist := {
    y0 <- nafill(dist, "locf")
    x0 <- nafill(replace(time, is.na(dist), NA), "locf")
    y1 <- nafill(dist, "nocb")
    x1 <- nafill(replace(time, is.na(dist), NA), "nocb")
    fifelse(is.na(dist), (y1 - y0) / (x1 - x0) * (time - x0) + y0, dist)
}]

временный код:

library(data.table)
set.seed(0L)
# df=data.frame(time=seq(7173,7195,1),dist=c(31091.33,NA,31100.00,31103.27,NA,NA,NA,NA,31124.98,NA,31132.81,NA,NA,NA,NA,31154.19,NA,31161.47,NA,NA,NA,NA,31182.97))
# DT=data.table(df)
nr <- 1e7
nNA <- nr/2
DT <- data.table(time=1:nr, dist=replace(rnorm(nr), sample(1:nr, nNA), NA_real_))

DT00 <- copy(DT)
DT01 <- copy(DT)
DT1 <- copy(DT)
DT20 <- copy(DT)
DT201 <- copy(DT)
DT202 <- copy(DT)
DT21 <- copy(DT)

mtd00 <- function() {
    DT00[, g := rleid(is.na(dist))]

    DT00[is.na(dist), dist := {
        i <- .I[c(1, .N)] + c(-1, 1)
        DT00[i[1]:i[2], approx(dist, n = .N)$y[-c(1, .N)]]
    }, by = g]
}

mtd01 <- function() {
    DT01[, g := rleid(is.na(dist))]

    DT01[is.na(dist), dist := {
        i <- .I[c(1, .N)] + c(-1, 1)
        DT01[i[1]:i[2], dist[1] + 1:(.N - 2)*(dist[.N] - dist[1])/(.N - 1)]
    }, by = g]
}

mtd1 <- function() {
    DT1[,dist_before := nafill(dist, "locf")]
    DT1[,dist_after := nafill(dist, "nocb")]
    DT1[, rle := rleid(dist)][,missings := max(.N +  1 , 2), by = rle][]
    DT1[is.na(dist), dist_before + .SD[,.I] *
            (dist_after - dist_before)/(missings), by = rle]
}


mtd20 <- function() {
    DT20[is.na(dist), {
        x0y0 <- DT20[!is.na(dist)][.SD, on=.(time), roll=Inf, .(time=x.time, dist=x.dist)]
        x1y1 <- DT20[!is.na(dist)][.SD, on=.(time), roll=-Inf, .(time=x.time, dist=x.dist)]
        (x1y1$dist - x0y0$dist) / (x1y1$time - x0y0$time) * (time - x0y0$time) + x0y0$dist
    }]
}

mtd201 <- function() {
    i <- DT201[, is.na(dist)]
    DT201[(i), {
        x0y0 <- DT201[(!i)][.SD, on=.(time), roll=Inf, .(time=x.time, dist=x.dist)]
        x1y1 <- DT201[(!i)][.SD, on=.(time), roll=-Inf, .(time=x.time, dist=x.dist)]
        (x1y1$dist - x0y0$dist) / (x1y1$time - x0y0$time) * (time - x0y0$time) + x0y0$dist
    }]
}

mtd202 <- function() {
    i <- DT201[is.na(dist), which=TRUE]
    DT201[i, {
        x0y0 <- DT201[-i][.SD, on=.(time), roll=Inf, .(time=x.time, dist=x.dist)]
        x1y1 <- DT201[-i][.SD, on=.(time), roll=-Inf, .(time=x.time, dist=x.dist)]
        (x1y1$dist - x0y0$dist) / (x1y1$time - x0y0$time) * (time - x0y0$time) + x0y0$dist
    }]
}


mtd21 <- function() {
    DT21[, {
        y0 <- nafill(dist, "locf")
        x0 <- nafill(replace(time, is.na(dist), NA), "locf")
        y1 <- nafill(dist, "nocb")
        x1 <- nafill(replace(time, is.na(dist), NA), "nocb")
        fifelse(is.na(dist), (y1 - y0) / (x1 - x0) * (time - x0) + y0, dist)
    }]
}

bench::mark(
    #mtd00(), mtd01(), 
    #mtd1(),
    mtd20(), mtd201(), mtd202(),
    mtd21(), check=FALSE)

тайминги:

# A tibble: 4 x 13
  expression      min   median `itr/sec` mem_alloc `gc/sec` n_itr  n_gc total_time result            memory            time    gc            
  <bch:expr> <bch:tm> <bch:tm>     <dbl> <bch:byt>    <dbl> <int> <dbl>   <bch:tm> <list>            <list>            <list>  <list>        
1 mtd20()       1.19s    1.19s     0.838    1.01GB    1.68      1     2      1.19s <dbl [5,000,000]> <df[,3] [292 x 3~ <bch:t~ <tibble [1 x ~
2 mtd201()      1.12s    1.12s     0.894  954.06MB    0.894     1     1      1.12s <dbl [5,000,000]> <df[,3] [341 x 3~ <bch:t~ <tibble [1 x ~
3 mtd202()      1.16s    1.16s     0.864  858.66MB    1.73      1     2      1.16s <dbl [5,000,000]> <df[,3] [392 x 3~ <bch:t~ <tibble [1 x ~
4 mtd21()    729.93ms 729.93ms     1.37   763.11MB    1.37      1     1   729.93ms <dbl [10,000,000~ <df[,3] [215 x 3~ <bch:t~ <tibble [1 x ~

отредактировать: обратиться к комментарию, используя is.na(dist)несколько раз:

set.seed(0L)
nr <- 1e7
nNA <- nr/2
DT <- data.table(time=1:nr, dist=replace(rnorm(nr), sample(1:nr, nNA), NA_real_))
DT20 <- copy(DT)
DT201 <- copy(DT)
DT202 <- copy(DT)

mtd20 <- function() {
    DT20[is.na(dist), dist := {
        x0y0 <- DT20[!is.na(dist)][.SD, on=.(time), roll=Inf, .(time=x.time, dist=x.dist)]
        x1y1 <- DT20[!is.na(dist)][.SD, on=.(time), roll=-Inf, .(time=x.time, dist=x.dist)]
        (x1y1$dist - x0y0$dist) / (x1y1$time - x0y0$time) * (time - x0y0$time) + x0y0$dist
    }]
}

mtd201 <- function() {
    i <- DT201[, is.na(dist)]
    DT201[(i), dist := {
        x0y0 <- DT201[(!i)][.SD, on=.(time), roll=Inf, .(time=x.time, dist=x.dist)]
        x1y1 <- DT201[(!i)][.SD, on=.(time), roll=-Inf, .(time=x.time, dist=x.dist)]
        (x1y1$dist - x0y0$dist) / (x1y1$time - x0y0$time) * (time - x0y0$time) + x0y0$dist
    }]
}

mtd202 <- function() {
    i <- DT201[is.na(dist), which=TRUE]
    DT201[i, dist := {
        x0y0 <- DT201[-i][.SD, on=.(time), roll=Inf, .(time=x.time, dist=x.dist)]
        x1y1 <- DT201[-i][.SD, on=.(time), roll=-Inf, .(time=x.time, dist=x.dist)]
        (x1y1$dist - x0y0$dist) / (x1y1$time - x0y0$time) * (time - x0y0$time) + x0y0$dist
    }]
}

тайминги:

# A tibble: 3 x 13
  expression      min   median `itr/sec` mem_alloc `gc/sec` n_itr  n_gc total_time result                    memory             time     gc               
  <bch:expr> <bch:tm> <bch:tm>     <dbl> <bch:byt>    <dbl> <int> <dbl>   <bch:tm> <list>                    <list>             <list>   <list>           
1 mtd20()      24.1ms   25.8ms      37.5    1.01GB    13.6     11     4      294ms <df[,2] [10,000,000 x 2]> <df[,3] [310 x 3]> <bch:tm> <tibble [15 x 3]>
2 mtd201()     24.8ms   25.6ms      38.2  954.07MB     8.19    14     3      366ms <df[,2] [10,000,000 x 2]> <df[,3] [398 x 3]> <bch:tm> <tibble [17 x 3]>
3 mtd202()       24ms   25.6ms      38.3   76.39MB     8.22    14     3      365ms <df[,2] [10,000,000 x 2]> <df[,3] [241 x 3]> <bch:tm> <tibble [17 x 3]>

Не вижу большой разницы в таймингах при уменьшении количества is.na(dist)звонков

chinsoon12
источник

1

is.na(dist)вычисляется 3 раза, его можно вычислить один раз и использовать повторно

jangorecki

Нелегко сравнивать время, когда есть смешанные единицы ( ms/ us)

jangorecki

Я не могу повторить результаты теста. DT_x <- copy(DT)вероятно, должен быть в начале каждого вызова функции. Обновление по ссылке происходит в вызовах функций.

Коул

@Cole спасибо, я всегда беспокоюсь, что копия повлияет на разницу во времени. поэтому я склонен оставлять это снаружи. что касается обновления по ссылке, 1) память уже выделена, 2) код не предполагает, что вычисленный столбец был предварительно вычислен, и не использует вычисляемый столбец, и 3) разбивка столбца происходит при каждом повторении и, следовательно, мы надеемся, что он имеет меньше влияние на сроки. для первого, вы могли бы хотеть времяbench::mark(copy(DT), copy(DT))

chinsoon12

1

Сроки во многом зависят от того, сколько есть АН. Первый вызов функции обновляется по ссылке и заменяет NA значениями. Все последующие звонки заменить нечего. Для 1e7примера, с copy(DT)с 27 мс, то mtd20()вызов принял 1.43s , используя копию и только 30 мс , если удалить копию из функции.

Коул

5

С помощью library(zoo)

DT[, dist := na.approx(dist)]

В качестве альтернативы, если вы предпочитаете придерживаться базовых функций R, а не использовать другой пакет, то вы можете сделать

DT[, dist := approx(.I, dist, .I)$y]

DWW
источник

5

Вот подход rcpp, который проходит через все один раз с дополнительным проходом для всех элементов NA.

Rcpp::sourceCpp(code = '
#include <Rcpp.h>
using namespace Rcpp;

// [[Rcpp::export]]
NumericVector rcpp_approx2D(IntegerVector x, NumericVector y) {
  double x_start = 0, y_start = 0, slope = 0;
  int count = 0;

  NumericVector y1 = clone(y); //added to not update-by-reference

  for(int i = 0; i < y1.size(); ++i){
    if (NumericVector::is_na(y1[i])){
      count++;
    } else {
      if (count != 0) {
        x_start = x[i-(count+1)];
        y_start = y1[i-(count+1)];
        slope = (y1[i] - y_start) / (x[i]- x_start);
        for (int j = 0; j < count; j++){
          y1[i-(count-j)] = y_start + slope * (x[i - (count - j)] - x_start);
        }
        count = 0;
      }
    }
  }
  return(y1);
}
')

Тогда в R:

DT[, rcpp_approx2D(time, dist)]

капуста
источник

data.table линейно интерполируя значения NA без групп

Ответы: