Растровая бизерианская корреляция - PullRequest
Новая
       44

Растровая бизерианская корреляция

0 голосов
/ 12 марта 2020

Есть ли способ выполнить бизериальную корреляцию или точечно-бизериальную корреляцию между тепловой картой и двоичным растром, используя QGIS, r или python, учитывая, что оба имеют одинаковый экстент, разрешение и CRS?

1 Ответ

0 голосов
/ 15 марта 2020

Я разобрался с методом. Используя R, я получаю значения как из дихотомического растра, так и из непрерывного. Таким образом, можно было рассчитать точечно-двоичную корреляцию через функцию biserial.cor . Тем не менее, поскольку я не эксперт в R, надеюсь, вы, ребята, можете сказать мне, если что-то выглядит не так. 

library(raster)
library(ltm)

pointb.correlation <- function(cnt_lyr_vector, dct_lyr_vector, output_dir)
{
  cat('\nPerforming point-biserial correlation ...\n\n')
  f_result <- list()

  for (i in 1:length(names(dct_lyr_vector))) {

    dct_lyr <- dct_lyr_vector[[i]]
    dct_lyr_name <- names(dct_lyr_vector[[i]])

    cat('Getting values from dichotomous layer', dct_lyr_name, '...\n')
    dct_lyr_values <- extract(dct_lyr, extent(dct_lyr))

    cat('Replacing missing values from', dct_lyr_name, '...\n\n')
    dct_lyr_values[is.na(dct_lyr_values)] <- 0

    cat('Getting correlation between', dct_lyr_name, 'and:', '\n\n')
    result <- list()

    for (j in 1:length(names(cnt_lyr_vector))) {

      start_time <-Sys.time()

      cnt_lyr <- cnt_lyr_vector[[j]]
      cnt_lyr_name <- names(cnt_lyr_vector[[j]])

      cat('->', cnt_lyr_name, '\n')

      cat('Getting values from category', cnt_lyr_name, '...\n')
      cnt_lyr_values <- extract(cnt_lyr, extent(cnt_lyr))

      cat('Replacing missing values from', cnt_lyr_name, '...\n')
      cnt_lyr_values[is.na(cnt_lyr_values)] <- 0

      cat('Doing the math, be patient :)', '\n')
      r_key <- paste(dct_lyr_name, cnt_lyr_name, sep = ".")
      result[[r_key]] <- biserial.cor(cnt_lyr_values, dct_lyr_values, use = c("complete.obs"), level = 1)

      end_time <-Sys.time()
      time_taken <- end_time - start_time
      cat('Time taken: ', time_taken,  '\n\n')
    }
    filename <- file.path(output_dir, dct_lyr_name)
    write.table(unlist(result), filename, row.names = TRUE, col.names=FALSE, sep = ",", eol = "\n")

    f_result <- append(f_result, result)
  }
  return(f_result)
}

exec <- function(dichotomous_lyr_dir, continuous_lyr_dir, output_dir)
{
  cnt_vector <- grep(".tif$", list.files(file.path(continuous_lyr_dir), all.files = F), ignore.case = TRUE, value = TRUE)
  cnt_vector_full_path <- file.path(continuous_lyr_dir, cnt_vector)
  cnt_stack <- raster::stack(cnt_vector_full_path)

  dct_vector <- grep(".tif$", list.files(file.path(dichotomous_lyr_dir), all.files = F), ignore.case = TRUE, value = TRUE)
  dct_vector_full_path <- file.path(dichotomous_lyr_dir, dct_vector)
  dct_stack <- raster::stack(dct_vector_full_path)

  corr_matrix <- pointb.correlation(cnt_stack, dct_stack, output_dir)

  return(corr_matrix)
}

dct_lyr_folder <- ''
cnt_folder <- ''
results_folder <- ''
matrix <- exec(dct_lyr_folder, cnt_folder, results_folder)
...