Замена NA в R с ближайшим значением

Ответ 1

Вот очень быстрый. Он использует findInterval, чтобы найти, какие две позиции следует учитывать для каждого NA в ваших исходных данных:

f1 <- function(dat) {
  N <- length(dat)
  na.pos <- which(is.na(dat))
  if (length(na.pos) %in% c(0, N)) {
    return(dat)
  }
  non.na.pos <- which(!is.na(dat))
  intervals  <- findInterval(na.pos, non.na.pos,
                             all.inside = TRUE)
  left.pos   <- non.na.pos[pmax(1, intervals)]
  right.pos  <- non.na.pos[pmin(N, intervals+1)]
  left.dist  <- na.pos - left.pos
  right.dist <- right.pos - na.pos

  dat[na.pos] <- ifelse(left.dist <= right.dist,
                        dat[left.pos], dat[right.pos])
  return(dat)
}

И здесь я проверяю его:

# sample data, suggested by @JeffAllen
dat <- as.integer(runif(50000, min=0, max=10))
dat[dat==0] <- NA

# computation times
system.time(r0 <- f0(dat))    # your function
# user  system elapsed 
# 5.52    0.00    5.52
system.time(r1 <- f1(dat))    # this function
# user  system elapsed 
# 0.01    0.00    0.03
identical(r0, r1)
# [1] TRUE

Ответ 2

Код ниже. Первоначальный вопрос не был полностью определен, я просил дать эти разъяснения:

• Гарантировано ли, что по крайней мере первая и/или последняя записи являются не-NA? нет [Нет]
• Что делать, если все записи в строке являются NA? [Оставить как есть]
• Вам небезразлично, как разделены связи, например, как обрабатывать средний NA в 1 3 NA NA NA 5 7? [Не заботьтесь/левая]
• У вас есть верхняя граница (S) на самом длинном смежном промежутке NA в строке? (Я думаю, что рекурсивное решение, если S мало, или решение dataframe с ifelse, если S велико, а число строк и столбцов велико.) [худший случай S может быть патологически большим, следовательно, рекурсия не следует использовать]

geoffjentry, вашим решением ваши узкие места будут серийным расчетом nearest.non.na.pos и последовательным назначением dat[na.pos] <- dat[non.na.pos[nearest.non.na.pos]] Для большого промежутка длины G все, что нам действительно нужно вычислить, это то, что первые (G/2, округлые) элементы заполняются слева, остальные справа. (Я мог бы опубликовать ответ, используя ifelse, но он будет похож.) Являются ли ваши критерии временем выполнения, эффективностью большого вывода, использованием временной памяти или четкости кода?

Возможные хитрости Coupla:

• нужно только вычислить N <- length(dat) один раз
• ускорение обычного случая: if (length(na.pos) == 0) пропустить строку, так как она не имеет NA
• if (length(na.pos) == length(dat)-1) (редкий) случай, когда есть только одна запись, отличная от NA, поэтому мы заполняем всю строку с ней

Контурное решение:

К сожалению, na.locf не работает на полном кадре данных, вы должны использовать sapply, row-wise:

na.fill_from_nn <- function(x) {
  row.na <- is.na(x)
  fillFromLeft <- na.locf(x, na.rm=FALSE) 
  fillFromRight <- na.locf(x, fromLast=TRUE, na.rm=FALSE)

  disagree <- rle(fillFromLeft!=fillFromRight)
  for (loc in (disagree)) { ...  resolve conflicts, row-wise }
}

sapply(dat, na.fill_from_nn)

В качестве альтернативы, так как, как вы говорите, смежные NA являются редкими, сделайте быстрый и немой ifelse, чтобы заполнить изолированные NAs слева. Это будет работать с data-frame wise = > делает обычный случай быстрым. Затем обрабатывайте все остальные случаи с помощью циклического цикла. (Это повлияет на тай-брейк на средних элементах в длинном промежутке НС, но вы говорите, что вам все равно.)

Ответ 3

Я не могу придумать очевидное простое решение, но, посмотрев на предложения (особенно smci предложение rle), я подошел со сложной функцией, которая оказывается более эффективной.

Это код, я объясню ниже:

# Your function
your.func = function(dat) {
  na.pos <- which(is.na(dat))
  if (length(na.pos) == length(dat)) {
    return(dat)
  }
  non.na.pos <- setdiff(seq_along(dat), na.pos)
  nearest.non.na.pos <- sapply(na.pos, function(x) which.min(abs(non.na.pos - x)))
  dat[na.pos] <- dat[non.na.pos[nearest.non.na.pos]]
  dat
}

# My function
my.func = function(dat) {
    nas=is.na(dat)
    if (!any(!nas)) return (dat)
    t=rle(nas)
    f=sapply(t$lengths[t$values],seq)
    a=unlist(f)
    b=unlist(lapply(f,rev))
    x=which(nas)
    l=length(dat)
    dat[nas]=ifelse(a>b,dat[ ifelse((x+b)>l,x-a,x+b) ],dat[ifelse((x-a)<1,x+b,x-a)])
    dat
}


# Test
n = 100000
test.vec = 1:n
set.seed(1)
test.vec[sample(test.vec,n/4)]=NA

system.time(t1<-my.func(test.vec))
system.time(t2<-your.func(test.vec)) # 10 times speed improvement on my machine

# Verify
any(t1!=t2)

Моя функция основана на rle. Я читаю комментарии выше, но мне кажется, что rle отлично работает для NA. Это проще всего объяснить небольшим примером.

Если я начинаю с вектора:

dat=c(1,2,3,4,NA,NA,NA,8,NA,10,11,12,NA,NA,NA,NA,NA,18)

Затем я получаю позиции всех NA:

x=c(5,6,7,8,13,14,15,16,17)

Затем для каждого "запуска" NA создайте последовательность от 1 до длины пробега:

a=c(1,2,3,1,1,2,3,4,5)

Затем я делаю это снова, но я меняю последовательность:

b=c(3,2,1,1,5,4,3,2,1)

Теперь я могу просто сравнить векторы a и b: если a <= b, оглянитесь назад и возьмите значение в x-a. Если a > b, то посмотрите вперед и возьмите значение в x + b. Остальное просто обрабатывает угловые случаи, когда у вас есть все NA или NA, которые работают в конце или в начале вектора.

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

Ответ 4

Здесь мой удар по нему. Я никогда не хотел бы видеть цикл for в R, но в случае редко-NA-вектора это выглядит действительно более эффективным (показатели производительности ниже). Суть кода ниже.

  #get the index of all NA values
  nas <- which(is.na(dat))

  #get the Boolean map of which are NAs, used later to determine which values can be used as a replacement, and which are just filled-in NA values
  namask <- is.na(dat)

  #calculate the maximum size of a run of NAs
  length <- getLengthNAs(dat);

  #the furthest away an NA value could be is half of the length of the maximum NA run
  windowSize <- ceiling(length/2)

  #loop through all NAs
  for (thisIndex in nas){
    #extract the neighborhood of this NA
    neighborhood <- dat[(thisIndex-windowSize):(thisIndex+windowSize)]
    #any already-filled-in values which were NA can be replaced with NAs
    neighborhood[namask[(thisIndex-windowSize):(thisIndex+windowSize)]] <- NA

    #the center of this neighborhood
    center <- windowSize + 1

    #compute the difference within this neighborhood to find the nearest non-NA value
    delta <- center - which(!is.na(neighborhood))

    #find the closest replacement
    replacement <- delta[abs(delta) == min(abs(delta))]
    #in case length > 1, just pick the first
    replacement <- replacement[1]

    #replace with the nearest non-NA value.
    dat[thisIndex] <- dat[(thisIndex - (replacement))]
  }

Мне понравился код, который вы предложили, но я заметил, что мы вычисляем дельту между каждым значением NA и каждым другим индексом, отличным от NA в матрице. Я думаю, что это был самый большой боевик. Вместо этого я просто извлекаю окрестности или окно минимального размера вокруг каждого NA и нахожу ближайшее не-NA значение в этом окне.

Таким образом, производительность масштабируется линейно по числу NA и размеру окна - где размер окна (потолок) составляет половину длины максимального пробега NA. Чтобы вычислить длину максимального пробега NA, вы можете использовать следующую функцию:

getLengthNAs <- function(dat){
  nas <- which(is.na(dat))
  spacing <- diff(nas)
  length <- 1;
  while (any(spacing == 1)){        
    length <- length + 1;
    spacing <- diff(which(spacing == 1))
  }
    length
}

Сравнение производительности

#create a test vector with 10% NAs and length 50,000.
dat <- as.integer(runif(50000, min=0, max=10))
dat[dat==0] <- NA

#the a() function is the code posted in the question
a <- function(dat){
  na.pos <- which(is.na(dat))
    if (length(na.pos) == length(dat)) {
        return(dat)
    }
    non.na.pos <- setdiff(seq_along(dat), na.pos)
    nearest.non.na.pos <- sapply(na.pos, function(x) {
        return(which.min(abs(non.na.pos - x)))
    })
    dat[na.pos] <- dat[non.na.pos[nearest.non.na.pos]]
    dat
}

#my code
b <- function(dat){
    #the same code posted above, but with some additional helper code to sanitize the input
    if(is.null(dat)){
      return(NULL);
    }

    if (all(is.na(dat))){
      stop("Can't impute NAs if there are no non-NA values.")
    }

    if (!any(is.na(dat))){
      return(dat);
    }

    #starts with an NA (or multiple), handle these
    if (is.na(dat[1])){
      firstNonNA <- which(!is.na(dat))[1]
      dat[1:(firstNonNA-1)] <- dat[firstNonNA]
    }

    #ends with an NA (or multiple), handle these
    if (is.na(dat[length(dat)])){
      lastNonNA <- which(!is.na(dat))
      lastNonNA <- lastNonNA[length(lastNonNA)]
      dat[(lastNonNA+1):length(dat)] <- dat[lastNonNA]
    }

    #get the index of all NA values
    nas <- which(is.na(dat))

    #get the Boolean map of which are NAs, used later to determine which values can be used as a replacement, and which are just filled-in NA values
    namask <- is.na(dat)

    #calculate the maximum size of a run of NAs
    length <- getLengthNAs(dat);

    #the furthest away an NA value could be is half of the length of the maximum NA run
    #if there a run at the beginning or end, then the nearest non-NA value could possibly be `length` away, so we need to keep the window large for that case.
    windowSize <- ceiling(length/2)

    #loop through all NAs
    for (thisIndex in nas){
      #extract the neighborhood of this NA
      neighborhood <- dat[(thisIndex-windowSize):(thisIndex+windowSize)]
      #any already-filled-in values which were NA can be replaced with NAs
      neighborhood[namask[(thisIndex-windowSize):(thisIndex+windowSize)]] <- NA

      #the center of this neighborhood
      center <- windowSize + 1

      #compute the difference within this neighborhood to find the nearest non-NA value
      delta <- center - which(!is.na(neighborhood))

      #find the closest replacement
      replacement <- delta[abs(delta) == min(abs(delta))]
      #in case length > 1, just pick the first
      replacement <- replacement[1]

      #replace with the nearest non-NA value.
      dat[thisIndex] <- dat[(thisIndex - (replacement))]
    }
    dat
}

#nograpes' answer on this question
c <- function(dat){
  nas=is.na(dat)
  if (!any(!nas)) return (dat)
  t=rle(nas)
  f=sapply(t$lengths[t$values],seq)
  a=unlist(f)
  b=unlist(lapply(f,rev))
  x=which(nas)
  l=length(dat)
  dat[nas]=ifelse(a>b,dat[ ifelse((x+b)>l,x-a,x+b) ],dat[ifelse((x-a)<1,x+b,x-a)])
  dat
}

#run 10 times each to get average performance.
sum <- 0; for (i in 1:10){ sum <- sum + system.time(a(dat))["elapsed"];}; cat ("A: ", sum/10)
A:  5.059
sum <- 0; for (i in 1:10){ sum <- sum + system.time(b(dat))["elapsed"];}; cat ("B: ", sum/10)
B:  0.126
sum <- 0; for (i in 1:10){ sum <- sum + system.time(c(dat))["elapsed"];}; cat ("C: ", sum/10)
C:  0.287

Таким образом, он выглядит как этот код (по крайней мере, в этих условиях), предлагает примерно 40-кратное ускорение от исходного кода, отправленного в вопросе, и ответ 2.2X по сравнению с @nograpes ниже (хотя я представляю себе rle решение, безусловно, будет быстрее в некоторых ситуациях, включая более богатый NA-вектор).

Ответ 5

Скорость примерно в 3-4 раза медленнее, чем у выбранного ответа. Мина довольно проста. Это редкий цикл while.

f2 <- function(x){

  # check if all are NA to skip loop
  if(!all(is.na(x))){

    # replace NA until they are gone
    while(anyNA(x)){

      # replace from the left
      x[is.na(x)] <- c(NA,x[1:(length(x)-1)])[is.na(x)]

      # replace from the right
      x[is.na(x)] <- c(x[-1],NA)[is.na(x)]
    }
  }

  # return original or fixed x
  x
}

Ответ 6

Мне нравятся все строгие решения. Хотя не напрямую, что было задано, я нашел это сообщение в поисках решения для заполнения значений NA интерполяцией. После рассмотрения этого сообщения я обнаружил na.fill на объекте зоопарка (вектор, фактор или матрицу):

z < -zoo (c (1,2,3,4,5,6, NA, NA, NA, 2,3,4,5,6, NA, NA, 4,6,7, NA ))

z1 < -na.fill(z, "extend" )

Обратите внимание на плавный переход по значениям NA

1.0 2.0 3.0 4.0 5.0 6.0 5.0 4.0 3.0 2.0 3.0 4.0 5.0 6.0 5.3 4.6 4.0 6.0 7.0 7.0

Возможно, это могло бы помочь