光栅双数相关

问题描述 投票:0回答:1

考虑到两者在程度,分辨率和CRS上都相同,是否可以通过使用QGIS,r或python在热图和二进制栅格之间执行双线性相关或点双线性相关?

python r correlation raster qgis
1个回答
0
投票

我想出了一种方法。通过使用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)
© www.soinside.com 2019 - 2024. All rights reserved.