#' Calculate the segregation index of NDiv
#'
#' The segregation index of NDiv assesses how much NDiv deviates from its random
#' expectation using a randomization approach. Values close to zero indicate a
#' random distribution of the traits that were used to calculate NDiv, values
#' larger than zero indicate a regular distribution (stem-wise mixing) and
#' values smaller than zero indicate a clustered configuration (patch-wise
#' mixing).
#'
#' @template apa_list_arg
#'
#' @param nsim Integer, how many simulation runs will be done.
#'
#' @param save_folder If specified, intermediate simulation results and
#' processing report will be stored in this folder (defaults to NULL).
#'
#' @param save_simulations If TRUE (default) and if \code{save_folder} is
#' specified, individual simulations will be stored (may consume large
#' storage capacities).
#'
#' @param overwrite If FALSE (default), \code{save_folder} needs to be empty. If
#' TRUE, all content in \code{save_folder} will be deleted before the new
#' simulations start.
#'
#' @param buffer_separate Should tree locations in the buffer and cores zones
#' be randomized separately? Defaults to TRUE.
#'
#' @param parallel If TRUE (default), parallel processing will be used.
#'
#' @param cl If \code{parallel} is TRUE, a Parallel Socket Cluster that was
#' created with the [parallel::makeCluster()]-function is specified here. If
#' NULL (default), \code{parallel::makeCluster(no_cores)} will be used.
#'
#' @param no_cores The number of cores that will be used for parallel processing
#' (default is 1). If NULL is specified, \code{parallel::detectCores() - 1}
#' will be used to automatically use one core less than available.
#'
#' @param iseed A seed that will be passed to [parallel::clusterSetRNGStream] to
#' make simulations reproducible when parallel processing is used (default is
#' 42).
#'
#' @details See Glatthorn (2021) for details.
#'
#' @template ref_glatthorn_2021
#'
#' @example examples/example_apa_seg.R
#'
#' @export
<- (,
nsim = 1000,
save_folder = ,
save_simulations = ,
overwrite = ,
buffer_separate = ,
= ,
cl = ,
no_cores = 1,
iseed = 42){
apa_config <- (, "apa_config")
apa_properties <- apa_config$apa_properties
(!(save_folder)){
(!(save_folder)){
(save_folder)
}
previous_runs <- (save_folder, full.names = )
(((save_folder)) != 0 & !overwrite){
("Directory `", save_folder, "` is not empty. Select `overwrite = TRUE`` to delete previous simulation runs.")
}
(previous_runs, recursive = )
}
sim_idx_nchar <- (nsim)
sim_idx <- (("%0", sim_idx_nchar, "d") , 1:nsim)
((apa_properties[1]) || !(("_ndiv", apa_properties)) ||
!(("_pdiv", apa_properties))){
("`ndiv` and `pdiv` are no an apa-properties in `apa_list`. Please use ",
"`apa_ndiv(apa_list, pdiv = TRUE)` to calculate ndiv before attempting to ",
"calculate the segregation of ndiv.")
}
apa_config <- (, "apa_config")
apa_properties <- apa_config$apa_properties
drop_tree_dat <- (x){
(!(x$tree_dat)){
x$tree_dat <- ()
}
x
}
((apa_config$dis_trait_column), "_ndiv")
keep_columns <- (apa_config$plot_id_column,
apa_config$tree_id_column,
apa_config$agg_class_column,
apa_config$subplot_id_column,
((apa_config$dis_trait_column), "_ndiv"),
((apa_config$dis_trait_column), "_pdiv"))
observed_output <-
lapply_deep(, .f = drop_columns_except, keep_columns = keep_columns,
.other = drop_tree_dat)
start_time <- ()
(){
((cl)){
((no_cores)){
no_cores <- ::() - 1
}
cl <- ::(no_cores)
}
# `loadNamespace("APAtree")` causes some Namespace problems that I completely fail to understand.
# `library("APAtree)` works fine.
(!(iseed)){
::(cl, iseed = iseed)
}
junk <- ::(cl, ("APAtree"))
sim_output <- ::(cl, sim_idx, nsim = nsim, fun = seg_randomize,
= , save_folder = save_folder, save_simulations = save_simulations,
buffer_separate = buffer_separate, keep_columns = keep_columns,
start_time = start_time)
::(cl)
}{
sim_output <- (sim_idx, nsim = nsim, FUN = seg_randomize,buffer_separate = buffer_separate,
= , save_folder = save_folder, save_simulations = save_simulations,
keep_columns = keep_columns,
start_time = start_time)
}
(!(save_folder)){
temp_files <- ("temp.txt", (save_folder, full.names = ), value = )
(temp_files)
}
apa_prop_columns <-
(((apa_config$dis_trait_column), "_ndiv"),
((apa_config$dis_trait_column), "_pdiv"))
sim_output <-
mapply_deep(((observed_output), sim_output),
combine_df_cols,
value_columns = apa_prop_columns)
append_to_column_name <- (x, column_names, ){
columns_in_x <- ((x), column_names)
new_column_names <- (columns_in_x, , sep = "_")
(x)[names(x) %in% column_names] <- new_column_names
x
}
sim_output <- lapply_deep(sim_output, .f = seg_orthogonalize)
sim_output <-
lapply_deep(sim_output, .f = append_to_column_name,
column_names = apa_prop_columns,
= "replications")
((save_folder)){
(sim_output, = (save_folder, "/sim_output.rda"))
}
seg_columns <- ("seg_", (apa_config$dis_trait_column), "_ndiv")
apa_list_seg <-
mapply_deep((, sim_output), .f = add_col,
= seg_columns)
(apa_list_seg) <- ("apa_list", "list")
apa_config_new <- apa_config
apa_config_new$apa_properties <- (apa_config_new$apa_properties, seg_columns)
(apa_list_seg, "apa_config") <- apa_config_new
validate_apa_list(apa_list_seg)
(!(save_folder)){
(apa_list_seg, = (save_folder, "/apa_list_seg.rda"))
}
apa_list_seg
}
#' @keywords internal
randomize_apa_list <-
(,
buffer_separate = ,
save_folder = ,
suffix = ,
= ){
apa_config <- (, "apa_config")
# Randomization is done over each each element in `plot_id_values`
((apa_config$plot_id_values) > 1){
apa_sub_list <- (X = apa_config$plot_id_values, FUN = ,
x = )
apa_sub_list_sim <- (apa_sub_list, randomize_apa_list, buffer_separate = buffer_separate,
save_folder = , suffix = , = )
apa_list_sim <- mapply_deep(apa_sub_list_sim, .f = (x){(what = , x)})
apa_config$randomized <-
(apa_list_sim, "apa_config") <- apa_config
}{
apa_list_sim <-
core_poly <- apa_list_sim$plot_dat[apa_config$core_column]]
buffer_poly <- apa_list_sim$plot_dat[apa_config$buffer_column]]
(!buffer_separate){
core_poly <- buffer_poly
}
donut_poly <- core_poly[NULL]
(poly_i (core_poly)){
donut_poly_i <- sf::st_difference(buffer_poly[poly_i],
core_poly[poly_i])
((donut_poly_i) != 0){
donut_poly <- (donut_poly, donut_poly_i)
}
}
donut_poly <- sf::st_sfc(donut_poly)
resample_poly <- (core_poly, donut_poly)
((sf::st_overlaps(resample_poly, sparse = ))){
("Some polygons for resampling are overlapping.")
}
intersect_mat <- sf::st_intersects(apa_list_sim$tree_dat, resample_poly, sparse = )
intersect_vec <- (intersect_mat)
intersect_vec!intersect_mat] <- 0
intersect_vec <- (intersect_vec, 1, )
resample_freq <- (intersect_vec)
(!"0" %in% (resample_freq)){
resample_freq <- ("0" = 0, resample_freq)
}
resample_points <-
(sf::st_geometry(apa_list_sim$tree_dat)[intersect_vec == 0],
sf::st_sample(resample_poly, size = resample_freq-1]))
resample_points_sorted <- resample_points[NULL]
resample_points_sorted[sort(intersect_vec, index.return = )$ix] <-
resample_points
sf::st_geometry(apa_list_sim$tree_dat) <- resample_points_sorted
apa_config$randomized <-
(apa_list_sim, "apa_config")$randomized <-
}
(){
apa_list_sim <- update_apa(apa_list_sim, = "all", no_polygon_update = )
}
(!(save_folder)){
object_name <- ("apa_list_sim_", suffix)
(object_name, apa_list_sim)
( = object_name,
= (save_folder, "/sim_", suffix, ".rda"))
}
validate_apa_list(apa_list_sim)
}
#' @keywords internal
seg_orthogonalize <-
(ndiv, pdiv = ) {
("seg_orthogonalize")
}
#' @keywords internal
seg_orthogonalize.default <-
(ndiv, pdiv){
ndiv_obs <- ndiv[1]
pdiv_obs <- pdiv[1]
((!(ndiv-1]) & !(pdiv-1])) < 3){
()
}
regression_data <-
(ndiv_sim = ndiv-1],
pdiv_sim = pdiv-1])
lm_1 <- ::(ndiv ~ pdiv, = regression_data)
pred <- ::(lm_1, newdata = (pdiv = pdiv_obs), se.fit = , interval = "prediction")
ndiv_mean_exp <- pred$fit[1, "fit"]
ndiv_sig_exp = (pred$se.fit^2 +pred$residual.scale^2)
(ndiv_obs - ndiv_mean_exp) / ndiv_sig_exp
}
#' @keywords internal
seg_orthogonalize.data.frame <- (ndiv){
ndiv_names <- (("_ndiv", (ndiv), value = ))
pdiv_names <- (("_pdiv", (ndiv), value = ))
(ndiv_i ndiv_names){
ndiv_i
trait_i <- ("_ndiv", "", ndiv_i)
pdiv_i <- (trait_i, "_pdiv")
seg_ndiv_i <- (seg_orthogonalize.default)(ndiv[ndiv_i]], ndiv[pdiv_i]])
ndiv[paste0("seg_", trait_i, "_ndiv")]] <- seg_ndiv_i
}
# ndiv <- ndiv[intersect(names(ndiv), c("species_ndiv", "species_pdiv"))]
# ndiv$seg_ndiv <- seg_ndiv
ndiv
}
#' @keywords internal
seg_randomize <-
(sim_i, nsim, , save_folder = , save_simulations = ,
buffer_separate = ,
keep_columns = keep_columns, start_time = ){
{
(!(save_folder)){
proc_file <- ((save_folder, "/", () , "_log.txt"),
= "at")
(proc_file, = "message")
}
("\n----------- SIMULATION ", sim_i, " -----------")
sim_output_i <- randomize_apa_list(, buffer_separate = buffer_separate,
save_folder = (save_simulations){save_folder}{},
suffix = sim_i)
(!(save_folder)){
( = "message")
(proc_file)
(sim_i, "\n", = (save_folder, "/", (), "_temp.txt"), = )
temp_files <- ("temp.txt", (save_folder, full.names = ), value = )
sim_proc <- (((temp_files, ), recursive = ))
proc_n <- (sim_proc)
time_passed <- ((), start_time)
time_total_est <- time_passed / proc_n * nsim
(
"Starting time of the simulations: ", (start_time), "\n",
proc_n, " of ", nsim , " simulations processed (",
(100 * proc_n / nsim, 2)," %).\n",
"Estimated time left: ", (time_total_est - time_passed, 2), " ",
(time_total_est, "units"), ".\n",
"Estimated finishing time: ",
(start_time + time_total_est), ".\n",
sep = "",
= (save_folder, "/proc.txt"))
}
drop_tree_dat <- (x){
(!(x"tree_dat"]])){
x$tree_dat <- ()
}
x
}
sim_output_i <-
lapply_deep(sim_output_i, .f = drop_columns_except,
keep_columns = keep_columns, .other = drop_tree_dat)
}
sim_output_i
}
