diff --git a/NAMESPACE b/NAMESPACE index af7799e9..d8d62b6a 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -51,6 +51,7 @@ export(network_listw.mc) export(nkde) export(nkde.mc) export(nkde_get_loo_values) +export(nkde_worker) export(quartic_kernel) export(simple_lines) export(simplify_network) diff --git a/R/adaptive_simultaneous_tnkde_bw.R b/R/adaptive_simultaneous_tnkde_bw.R index 6ae9ba1d..793962c7 100644 --- a/R/adaptive_simultaneous_tnkde_bw.R +++ b/R/adaptive_simultaneous_tnkde_bw.R @@ -152,7 +152,7 @@ adaptive_bw_tnkde <- function(grid, events_loc, events, lines, tol, digits, sparse, verbose){ ##step 1 split the datas ! - selections <- split_by_grid_abw(grid, events_loc, lines, max(trim_bw_net), tol, digits) + selections <- split_by_grid_abw.mc(grid, events_loc, lines, max(trim_bw_net), tol, digits) ## step 2 calculating the temp NKDE values if(verbose){ @@ -200,6 +200,11 @@ adaptive_bw_tnkde <- function(grid, events_loc, events, lines, tot_arr <- do.call(abind::abind, lapply(dfs, function(x){x[[1]]})) all_wids <- do.call(c, lapply(dfs, function(x){x[[2]]})) + ## petite galipette pour eviter des valeurs incroyablement basses ! + tot_arr <- ifelse(tot_arr < 2.225074e-200, 0, tot_arr) + above_0 <- min(c(tot_arr)[c(tot_arr) > 0]) + tot_arr <- ifelse(tot_arr == 0, above_0, tot_arr) + # tot_df <- do.call(rbind,dfs) # tot_df <- data.frame(tot_df[order(tot_df[,1]),]) # names(tot_df) <- c("goid", "k") @@ -354,6 +359,11 @@ adaptive_bw_tnkde.mc <- function(grid, events_loc, events, lines, tot_arr <- do.call(abind::abind, lapply(dfs, function(x){x[[1]]})) all_wids <- do.call(c, lapply(dfs, function(x){x[[2]]})) + ## petite galipette pour eviter des valeurs incroyablement basses ! + tot_arr <- ifelse(tot_arr < 2.225074e-200, 0, tot_arr) + above_0 <- min(c(tot_arr)[c(tot_arr) > 0]) + tot_arr <- ifelse(tot_arr == 0, above_0, tot_arr) + # tot_df <- do.call(rbind,dfs) # tot_df <- data.frame(tot_df[order(tot_df[,1]),]) # names(tot_df) <- c("goid", "k") diff --git a/R/nkde_execution_functions_sf.R b/R/nkde_execution_functions_sf.R index 47244ac1..3e9da394 100644 --- a/R/nkde_execution_functions_sf.R +++ b/R/nkde_execution_functions_sf.R @@ -858,7 +858,7 @@ adaptive_bw.mc <- function(grid,events,lines,bw,trim_bw,method,kernel_name,max_d rep(x, nrow(sel$events)) }) } - invisible(capture.output(values <- nkde_worker(lines = sel$lines, + invisible(capture.output(values <- spNetwork::nkde_worker(lines = sel$lines, events = sel$events, samples = sel$samples, kernel_name = kernel_name, @@ -900,7 +900,7 @@ adaptive_bw.mc <- function(grid,events,lines,bw,trim_bw,method,kernel_name,max_d }) } - values <- nkde_worker(lines = sel$lines, + values <- spNetwork::nkde_worker(lines = sel$lines, events = sel$events, samples = sel$samples, kernel_name = kernel_name, @@ -996,6 +996,7 @@ adaptive_bw.mc <- function(grid,events,lines,bw,trim_bw,method,kernel_name,max_d #' @importFrom igraph adjacent_vertices get.edge.ids #' @return A numeric vector with the nkde values #' @keywords internal +#' @export #' @examples #' #This is an internal function, no example provided nkde_worker <- function(lines, events, samples, kernel_name, bw, bws, method, div, digits, tol, sparse, max_depth, verbose = FALSE){ diff --git a/R/temporal_nkde_sf.R b/R/temporal_nkde_sf.R index ce21f296..4e7a3a9d 100644 --- a/R/temporal_nkde_sf.R +++ b/R/temporal_nkde_sf.R @@ -483,11 +483,16 @@ tnkde <- function(lines, events, time_field, w, samples_loc, samples_time, kerne }else{ if(adaptive_separate == TRUE){ ## we want to use an adaptive bw in the network space - bws_net <- adaptive_bw(grid, events_loc, lines, bw_net, trim_bw_net, method, + bws_net_all <- adaptive_bw.mc(grid = grid, + events = events_loc, + lines = lines, + bw = bw_net, + trim_bw = trim_bw_net, + method, kernel_name, max_depth, tol, digits, sparse, verbose) ## and in the time space bws_time <- adaptive_bw_1d(events$time, w, bw_time, kernel_name) - bws_net <- bws_net[events$goid] + bws_net <- bws_net_all[events$goid] }else{ interaction_bws <- adaptive_bw_tnkde(grid, events_loc, events, lines, bw_net, bw_time ,trim_bw_net, trim_bw_time, diff --git a/docs/articles/Isochrones.html b/docs/articles/Isochrones.html index 403c46e0..c93ac9f4 100644 --- a/docs/articles/Isochrones.html +++ b/docs/articles/Isochrones.html @@ -106,7 +106,7 @@

Calculating isochrones

Jeremy Gelb

-

2023-11-01

+

2024-01-01

Source: vignettes/Isochrones.Rmd diff --git a/docs/articles/KNetworkFunctions.html b/docs/articles/KNetworkFunctions.html index 44a369c6..4fa7bc33 100644 --- a/docs/articles/KNetworkFunctions.html +++ b/docs/articles/KNetworkFunctions.html @@ -106,7 +106,7 @@

Network k Functions

Jeremy Gelb

-

2023-11-01

+

2024-01-01

Source: vignettes/KNetworkFunctions.Rmd diff --git a/docs/articles/KNetworkFunctions_files/figure-html/unnamed-chunk-10-1.png b/docs/articles/KNetworkFunctions_files/figure-html/unnamed-chunk-10-1.png index 2d08f868..b510321c 100644 Binary files a/docs/articles/KNetworkFunctions_files/figure-html/unnamed-chunk-10-1.png and b/docs/articles/KNetworkFunctions_files/figure-html/unnamed-chunk-10-1.png differ diff --git a/docs/articles/KNetworkFunctions_files/figure-html/unnamed-chunk-12-1.png b/docs/articles/KNetworkFunctions_files/figure-html/unnamed-chunk-12-1.png index d64abd4d..4b420123 100644 Binary files a/docs/articles/KNetworkFunctions_files/figure-html/unnamed-chunk-12-1.png and b/docs/articles/KNetworkFunctions_files/figure-html/unnamed-chunk-12-1.png differ diff --git a/docs/articles/KNetworkFunctions_files/figure-html/unnamed-chunk-14-1.png b/docs/articles/KNetworkFunctions_files/figure-html/unnamed-chunk-14-1.png index 9f1b71f0..a594fc23 100644 Binary files a/docs/articles/KNetworkFunctions_files/figure-html/unnamed-chunk-14-1.png and b/docs/articles/KNetworkFunctions_files/figure-html/unnamed-chunk-14-1.png differ diff --git a/docs/articles/KNetworkFunctions_files/figure-html/unnamed-chunk-4-1.png b/docs/articles/KNetworkFunctions_files/figure-html/unnamed-chunk-4-1.png index 0a7aca7c..0533a932 100644 Binary files a/docs/articles/KNetworkFunctions_files/figure-html/unnamed-chunk-4-1.png and b/docs/articles/KNetworkFunctions_files/figure-html/unnamed-chunk-4-1.png differ diff --git a/docs/articles/KNetworkFunctions_files/figure-html/unnamed-chunk-6-1.png b/docs/articles/KNetworkFunctions_files/figure-html/unnamed-chunk-6-1.png index 79f056a0..3a0171c2 100644 Binary files a/docs/articles/KNetworkFunctions_files/figure-html/unnamed-chunk-6-1.png and b/docs/articles/KNetworkFunctions_files/figure-html/unnamed-chunk-6-1.png differ diff --git a/docs/articles/KNetworkFunctions_files/figure-html/unnamed-chunk-8-1.png b/docs/articles/KNetworkFunctions_files/figure-html/unnamed-chunk-8-1.png index a34b3a4d..780647c8 100644 Binary files a/docs/articles/KNetworkFunctions_files/figure-html/unnamed-chunk-8-1.png and b/docs/articles/KNetworkFunctions_files/figure-html/unnamed-chunk-8-1.png differ diff --git a/docs/articles/NKDE.html b/docs/articles/NKDE.html index 6012a9e3..397a9258 100644 --- a/docs/articles/NKDE.html +++ b/docs/articles/NKDE.html @@ -106,7 +106,7 @@

Network Kernel Density Estimate

Jeremy Gelb

-

2023-11-01

+

2024-01-01

Source: vignettes/NKDE.Rmd @@ -846,24 +846,24 @@

Bandwidth selection
 bws_selection_cv <- bw_cv_likelihood_calc(
-  bws = seq(200,700,50),
+  bws = seq(50,700,50),
   lines = mtl_network, events = bike_accidents,
   w = rep(1,nrow(bike_accidents)),
   kernel_name = "quartic", method = "discontinuous",
   diggle_correction = FALSE, study_area = NULL,
   max_depth = 8,
-  digits=2, tol=0.1, agg=5, 
+  digits=2, tol=0.1, agg=5,
   sparse=TRUE, grid_shape=c(1,1),
   verbose=FALSE, check=TRUE)
 
 bws_selection_cvl <- bw_cvl_calc(
-  bws = seq(50,700, 50),
+  bws = seq(50,700,50),
   lines = mtl_network, events = bike_accidents,
   w = rep(1,nrow(bike_accidents)),
   kernel_name = "quartic", method = "discontinuous",
   diggle_correction = FALSE, study_area = NULL,
   max_depth = 8,
-  digits=2, tol=0.1, agg=5, 
+  digits=2, tol=0.1, agg=5,
   sparse=TRUE, grid_shape=c(1,1),
   verbose=FALSE, check=TRUE)

It is also possible to do bandwidth selection with adaptive @@ -872,22 +872,52 @@

Bandwidth selection
-bws_selection_cv_adpt <- bw_cv_likelihood_calc(
-  bws = seq(200,700,50),
-  lines = mtl_network, events = bike_accidents,
+bws_selection_cv_adpt_dis <- bw_cv_likelihood_calc(
+  bws = seq(50,700,50),
+  trim_bws = seq(50,700,50)*2,
+  lines = mtl_network,
+  events = bike_accidents,
   w = rep(1,nrow(bike_accidents)),
   adaptive = TRUE,
-  kernel_name = "quartic", method = "discontinuous",
-  diggle_correction = FALSE, study_area = NULL,
+  kernel_name = "quartic",
+  method = "discontinuous",
+  diggle_correction = FALSE,
+  study_area = NULL,
   max_depth = 8,
-  digits=2, tol=0.1, agg=5, 
-  sparse=TRUE, grid_shape=c(1,1),
-  verbose=FALSE, check=TRUE)
+  digits=2,
+  tol=0.1,
+  agg=5,
+  sparse=TRUE,
+  grid_shape=c(1,1),
+  verbose=TRUE,
+  check=TRUE)
+
+bws_selection_cv_adpt_cont <- bw_cv_likelihood_calc(
+  bws = seq(50,700,50),
+  trim_bws = seq(50,700,50)*2,
+  lines = mtl_network,
+  events = bike_accidents,
+  w = rep(1,nrow(bike_accidents)),
+  adaptive = TRUE,
+  kernel_name = "quartic",
+  method = "continuous",
+  diggle_correction = FALSE,
+  study_area = NULL,
+  max_depth = 8,
+  digits=2,
+  tol=0.1,
+  agg=5,
+  sparse=TRUE,
+  grid_shape=c(1,1),
+  verbose=TRUE,
+  check=TRUE)
 
 cv_values <- data.frame(
   "bw" = bws_selection_cv$bw,
   "cv_likelihood" = bws_selection_cv$cv_scores,
-  "cvl_crit" = bws_selection_cvl$cvl_scores
+  "cvl_crit" = bws_selection_cvl$cvl_scores,
+  "cv_likelihood_adpt_dis" = bws_selection_cv_adpt_dis$cv_scores,
+  "cv_likelihood_adpt_cont" = bws_selection_cv_adpt_cont$cv_scores
 )
@@ -901,7 +931,10 @@

Bandwidth selection -cv_likelihood_adpt +cv_likelihood_adpt_dis + +

@@ -918,6 +951,9 @@

Bandwidth selection -616.01 +

+ + + + + + + + + + + + +
+cv_likelihood_adpt_cont
+-391.58 +
@@ -932,6 +968,9 @@

Bandwidth selection -361.06

+-211.81 +
@@ -946,6 +985,9 @@

Bandwidth selection -158.79

+-125.84 +
@@ -960,6 +1002,9 @@

Bandwidth selection -86.49

+-69.96 +
@@ -974,6 +1019,9 @@

Bandwidth selection -48.56

+-32.20 +
@@ -988,6 +1036,9 @@

Bandwidth selection -28.58

+-22.27 +
@@ -1002,6 +1053,9 @@

Bandwidth selection -24.57

+-20.38 +
@@ -1016,6 +1070,9 @@

Bandwidth selection -22.58

+-18.51 +
@@ -1030,6 +1087,9 @@

Bandwidth selection -14.69

+-16.64 +
@@ -1044,6 +1104,9 @@

Bandwidth selection -14.69

+-16.74 +
@@ -1058,6 +1121,9 @@

Bandwidth selection -14.73

+-14.85 +
@@ -1072,6 +1138,9 @@

Bandwidth selection -14.78

+-14.93 +
@@ -1086,6 +1155,9 @@

Bandwidth selection -14.83

+-15.02 +
@@ -1100,6 +1172,9 @@

Bandwidth selection -14.89

+-15.10 +
@@ -1107,7 +1182,8 @@

Bandwidth selection

Complementary functions diff --git a/docs/articles/NetworkBuilding.html b/docs/articles/NetworkBuilding.html index fe0ca210..aba2617f 100644 --- a/docs/articles/NetworkBuilding.html +++ b/docs/articles/NetworkBuilding.html @@ -106,7 +106,7 @@

Building graphs

Jeremy Gelb

-

2023-11-01

+

2024-01-01

Source: vignettes/NetworkBuilding.Rmd @@ -286,13 +286,13 @@

Calculating for e ## Pearson's product-moment correlation ## ## data: bike_accidents$density and bike_accidents$btws -## t = 1.4765, df = 345, p-value = 0.1407 +## t = NA, df = 345, p-value = NA ## alternative hypothesis: true correlation is not equal to 0 ## 95 percent confidence interval: -## -0.02625894 0.18299854 +## NA NA ## sample estimates: -## cor -## 0.0792427 +## cor +## NA

We can see that there is no correlation between the two variables. There is no association between the degree of centrality of an accident location in the network and the density of accidents in a radius of 300m diff --git a/docs/articles/NetworkBuilding_files/figure-html/unnamed-chunk-11-1.png b/docs/articles/NetworkBuilding_files/figure-html/unnamed-chunk-11-1.png index 855e0b0e..fa7d4b47 100644 Binary files a/docs/articles/NetworkBuilding_files/figure-html/unnamed-chunk-11-1.png and b/docs/articles/NetworkBuilding_files/figure-html/unnamed-chunk-11-1.png differ diff --git a/docs/articles/NetworkBuilding_files/figure-html/unnamed-chunk-11-2.png b/docs/articles/NetworkBuilding_files/figure-html/unnamed-chunk-11-2.png index f6ae5df7..f24a8eb5 100644 Binary files a/docs/articles/NetworkBuilding_files/figure-html/unnamed-chunk-11-2.png and b/docs/articles/NetworkBuilding_files/figure-html/unnamed-chunk-11-2.png differ diff --git a/docs/articles/SpatialWeightMatrices.html b/docs/articles/SpatialWeightMatrices.html index a0e77688..2af9e629 100644 --- a/docs/articles/SpatialWeightMatrices.html +++ b/docs/articles/SpatialWeightMatrices.html @@ -106,7 +106,7 @@

Spatial Weight Matrices

Jeremy Gelb

-

2023-11-01

+

2024-01-01

Source: vignettes/SpatialWeightMatrices.Rmd diff --git a/docs/articles/TNKDE.html b/docs/articles/TNKDE.html index 79a45f4f..085e1997 100644 --- a/docs/articles/TNKDE.html +++ b/docs/articles/TNKDE.html @@ -106,7 +106,7 @@

Temporal Network Kernel Density Estimate

Jeremy Gelb

-

2023-11-01

+

2024-01-01

Source: vignettes/TNKDE.Rmd @@ -521,17 +521,15 @@

Adaptive bandwidth
-cv_scores_adpt <- bw_tnkde_cv_likelihood_calc(
-  bw_net_range = c(200,1100),
-  bw_net_step = 100,
-  bw_time_range = c(10,70),
-  bw_time_step = 10,
+cv_scores_adpt_dis <- bw_tnkde_cv_likelihood_calc(
+  bws_net = seq(200,1100,100),
+  bws_time = seq(10,70,10),
   lines = mtl_network,
   events = bike_accidents,
   time_field = "Time",
   adaptive = TRUE,
-  trim_net_bws = seq(200,1100,100),
-  trim_time_bws = seq(10,70,10),
+  trim_net_bws = seq(200,1100,100)*2,
+  trim_time_bws = seq(10,70,10)*2,
   w = rep(1, nrow(bike_accidents)),
   kernel_name = "quartic",
   method = "discontinuous",
@@ -545,9 +543,33 @@ 

Adaptive bandwidth=c(1,1), sub_sample=1, verbose = FALSE, + check = TRUE) + +cv_scores_adpt_cont <- bw_tnkde_cv_likelihood_calc.mc( + bws_net = seq(200,1100,100), + bws_time = seq(10,70,10), + lines = mtl_network, + events = bike_accidents, + time_field = "Time", + adaptive = TRUE, + trim_net_bws = seq(200,1100,100)*2, + trim_time_bws = seq(10,70,10)*2, + w = rep(1, nrow(bike_accidents)), + kernel_name = "quartic", + method = "continuous", + diggle_correction = FALSE, + study_area = NULL, + max_depth = 10, + digits = 2, + tol = 0.1, + agg = 15, + sparse=TRUE, + grid_shape=c(2,2), + sub_sample=1, + verbose = FALSE, check = TRUE)

-knitr::kable(cv_scores_adpt)
+knitr::kable(cv_scores_adpt_cont) @@ -572,108 +594,108 @@

Adaptive bandwidth

- - - - - - - + + + + + + + - - - - - - - + + + + + + + - - - - - - - + + + + + + + - - - - - - - + + + + + + + - - - - - - - + + + + + + + - - - - - - - + + + + + + + - - - - - - - + + + + + + + - - - - - - - + + + + + + + - - - - - - - + + + + + + + - - - - - - - + + + + + + +
200-610.29113-508.88811-469.66223-431.20407-400.65001-362.10526-317.80491-338.97156-292.67876-256.54268-224.32854-196.12899-165.89033-151.78938
300-508.22575-386.54432-329.30490-264.55071-230.06896-189.65467-169.32609-272.11248-189.42585-133.06588-104.92356-78.75489-56.66857-56.73603
400-434.90614-272.46528-215.60407-157.08051-112.41486-94.16037-81.97355-211.52817-108.49105-74.37144-52.34417-38.33608-36.44765-32.55607
500-327.16538-172.99572-114.17486-65.74005-51.37540-35.25870-33.20466-148.83050-68.21414-40.13705-34.23911-30.33735-26.47279-26.61781
600-227.77256-112.01773-73.51592-47.26551-41.07479-28.98670-26.99484-98.28102-56.15369-38.22042-32.35431-26.49237-24.63671-24.78364
700-164.78403-81.48058-53.20828-41.04888-32.91215-26.89284-26.92499-69.96943-48.15285-36.35144-28.46047-26.64666-24.78786-24.92029
800-121.77890-57.09647-40.93748-34.84482-30.84118-24.86597-24.92302-67.94635-38.15241-30.41514-28.55641-24.73862-24.89476-25.03609
900-79.28756-44.74426-34.73286-30.76957-26.82098-24.88398-24.95458-57.87577-34.18940-30.45316-26.63735-24.84286-25.00491-25.15673
1000-70.76108-40.59921-32.64436-28.74543-26.82357-24.91819-25.00168-55.91861-34.27345-30.54280-26.74034-24.94651-25.11732-25.27386
1100-66.59316-38.50578-32.61891-26.74280-26.85349-24.96216-25.05648-53.99006-32.35650-30.63393-30.90249-25.05491-25.22689-25.38583
-

The best score is reached with a network and temporal bandwidths of -800m and 60 days.

+

For a discontinuous NKDE The best score is reached with a network and +temporal bandwidths of 800m and 60 days.

 tnkde_densities <- tnkde(lines = mtl_network,
                    events = bike_accidents,
diff --git a/docs/articles/images/animated_map.gif b/docs/articles/images/animated_map.gif
index a85183c1..81825931 100644
Binary files a/docs/articles/images/animated_map.gif and b/docs/articles/images/animated_map.gif differ
diff --git a/docs/articles/web_vignettes/AdaptiveBW.html b/docs/articles/web_vignettes/AdaptiveBW.html
index 64b289f5..85e08fc0 100644
--- a/docs/articles/web_vignettes/AdaptiveBW.html
+++ b/docs/articles/web_vignettes/AdaptiveBW.html
@@ -106,7 +106,7 @@ 

K-nearest neighbour adaptive bandwidth

Jeremy Gelb

-

2023-11-01

+

2024-01-01

Source: vignettes/web_vignettes/AdaptiveBW.Rmd @@ -283,8 +283,9 @@

K nearest neighbours a #> Animation saved to E:\R\Packages\spNetwork\vignettes\web_vignettes\images\animated_densities.gif

 knitr::include_graphics("images/animated_densities.gif")
-

-With no suprise, the larger the bandwidth, the smoother the results

+

+

With no suprise, the larger the bandwidth, the smoother the +results

Finding the optimal bandwidth

@@ -299,7 +300,7 @@

Finding the optimal bandwidthnearest_nb_mat <- knn_dists$distances[,5:25] colnames(nearest_nb_mat) <- 5:25 -bw_scores <- bw_cv_likelihood_calc(lines = mtl_network, +bw_scores_dis <- bw_cv_likelihood_calc(lines = mtl_network, events = bike_accidents_agg, w = bike_accidents_agg$weight, kernel_name = 'quartic', @@ -307,31 +308,82 @@

Finding the optimal bandwidth adaptive = TRUE, mat_bws = nearest_nb_mat, digits = 1, - tol = 0.01, - max_depth = 10, agg = NULL, sparse = TRUE) -#> [1] "checking inputs ..." -#> [1] "prior data preparation ..." -#> [1] "Calculating the correction factor if required" -#> [1] "Splittig the dataset within the grid ..." -#> [1] "start calculating the CV values ..." -#> - | - | | 0% - | - |======================================================================| 100% - -df <- data.frame( - k = 5:25, - score = bw_scores[,2] -) -ggplot(df) + - geom_line(aes(x = k, y = score)) + - geom_point(aes(x = k, y = score), color = 'red') + - theme_bw()

+ tol = 0.01, + grid_shape = c(2,2), + max_depth = 8, agg = NULL, sparse = TRUE) +#> [1] "checking inputs ..." +#> [1] "prior data preparation ..." +#> [1] "Calculating the correction factor if required" +#> [1] "Splittig the dataset within the grid ..." +#> [1] "start calculating the CV values ..." +#> + | + | | 0%[1] 1 +#> + | + |================== | 25%[1] 2 +#> + | + |=================================== | 50%[1] 3 +#> + | + |==================================================== | 75%[1] 4 +#> + | + |======================================================================| 100% + +bw_scores_cont <- bw_cv_likelihood_calc(lines = mtl_network, + events = bike_accidents_agg, + w = bike_accidents_agg$weight, + kernel_name = 'quartic', + method = 'continuous', + adaptive = TRUE, + mat_bws = nearest_nb_mat, + digits = 1, + tol = 0.01, + grid_shape = c(2,2), + max_depth = 8, agg = NULL, sparse = TRUE) +#> [1] "checking inputs ..." +#> [1] "prior data preparation ..." +#> [1] "Calculating the correction factor if required" +#> [1] "Splittig the dataset within the grid ..." +#> [1] "start calculating the CV values ..." +#> + | + | | 0%[1] 1 +#> + | + |================== | 25%[1] 2 +#> + | + |=================================== | 50%[1] 3 +#> + | + |==================================================== | 75%[1] 4 +#> + | + |======================================================================| 100% + +df <- data.frame( + k = 5:25, + score_dis = bw_scores_dis[,2], + score_cont = bw_scores_cont[,2] +) + +ggplot(df) + + geom_line(aes(x = k, y = score_dis)) + + geom_point(aes(x = k, y = score_dis, color = 'discontinuous')) + + geom_line(aes(x = k, y = score_cont)) + + geom_point(aes(x = k, y = score_cont, color = 'continuous')) + + scale_color_manual(values = c("discontinuous" = "red", + "continuous" = "blue")) + + labs(color = '') + + theme_bw()

Based on these results, we decide to use the 8th nearest neighbour as -an optimal bandwidth. It clearly minimizes the likelihood cross -validation criterion.

+an optimal bandwidth. It clearly maximizes the likelihood cross +validation criterion for both the continuous and the discontinuous +NKDE.

@@ -481,7 +533,7 @@

Finding the optimal bandwidth} -values <- bw_tnkde_cv_likelihood_calc(lines = mtl_network, +scores_dis <- bw_tnkde_cv_likelihood_calc(lines = mtl_network, events = bike_accidents, time_field = "int_time", w = rep(1,nrow(bike_accidents)), @@ -505,18 +557,49 @@

Finding the optimal bandwidth=FALSE, check=TRUE) +scores_cont <- bw_tnkde_cv_likelihood_calc(lines = mtl_network, + events = bike_accidents, + time_field = "int_time", + w = rep(1,nrow(bike_accidents)), + kernel_name = 'quartic', + method = 'continuous', + arr_bws_net = network_bws, + arr_bws_time = time_bws, + diggle_correction = FALSE, + study_area = NULL, + adaptive = TRUE, + trim_net_bws = NULL, + trim_time_bws = NULL, + max_depth = 10, + digits=5, + tol=0.1, + agg=NULL, + sparse=TRUE, + zero_strat = "min_double", + grid_shape=c(1,1), + sub_sample=1, + verbose=FALSE, + check=TRUE) + df <- data.frame( k = 5:25, - score = values[,1] + score_dis = scores_dis[,1], + score_cont = scores_cont[,1] ) + ggplot(df) + - geom_line(aes(x = k, y = score)) + - geom_point(aes(x = k, y = score), color = 'red') + + geom_line(aes(x = k, y = score_dis)) + + geom_point(aes(x = k, y = score_dis, color = 'discontinuous')) + + geom_line(aes(x = k, y = score_cont)) + + geom_point(aes(x = k, y = score_cont, color = 'continuous')) + + scale_color_manual(values = c("discontinuous" = "red", + "continuous" = "blue")) + + labs(color = '') + theme_bw()

-

By checking the table above, we decide to select the 10th nearest -neighbour because it clearly corresponds to the major elbow in the above -figure. Also, the 21th neighbour gives the highest score.

+

By checking the table above, we could decide to select the 10th +nearest neighbour because it clearly corresponds to the major elbow in +the above figure. Also, the 21th neighbour gives the highest score.

Abramson, Ian S. 1982. “On Bandwidth Variation in Kernel diff --git a/docs/articles/web_vignettes/AdaptiveBW_files/figure-html/unnamed-chunk-14-1.png b/docs/articles/web_vignettes/AdaptiveBW_files/figure-html/unnamed-chunk-14-1.png index 9e61605e..0838bffd 100644 Binary files a/docs/articles/web_vignettes/AdaptiveBW_files/figure-html/unnamed-chunk-14-1.png and b/docs/articles/web_vignettes/AdaptiveBW_files/figure-html/unnamed-chunk-14-1.png differ diff --git a/docs/articles/web_vignettes/AdaptiveBW_files/figure-html/unnamed-chunk-8-1.png b/docs/articles/web_vignettes/AdaptiveBW_files/figure-html/unnamed-chunk-8-1.png index b8973d14..70bb0a40 100644 Binary files a/docs/articles/web_vignettes/AdaptiveBW_files/figure-html/unnamed-chunk-8-1.png and b/docs/articles/web_vignettes/AdaptiveBW_files/figure-html/unnamed-chunk-8-1.png differ diff --git a/docs/articles/web_vignettes/NKDEdetailed.html b/docs/articles/web_vignettes/NKDEdetailed.html index 96fae926..58f32fce 100644 --- a/docs/articles/web_vignettes/NKDEdetailed.html +++ b/docs/articles/web_vignettes/NKDEdetailed.html @@ -107,7 +107,7 @@

Details about NKDE

Jeremy Gelb

-

2023-11-01

+

2024-01-01

Source: vignettes/web_vignettes/NKDEdetailed.Rmd @@ -191,8 +191,8 @@

The simple method
 d3_plot_situation(all_lines, event, samples_pts, simple_kernel, scales = c(1,1,3))

-
-

(The function used to display the 3d plot is available in the source +

+

(The function used to display the 3d plot is available in the source code of the vignette)

As one can see, NKDE is only determined by the distance between the sampling points and the event. This is a problem, because at @@ -247,8 +247,8 @@

The discontinuous method verbose = FALSE)

 d3_plot_situation(all_lines, event, samples_pts, discontinuous_kernel, scales = c(1,1,3))
-
-

As one can see, the values of the NKDE are split at intersections to +

+

As one can see, the values of the NKDE are split at intersections to avoid the multiplication of the mass observed in the simple version. However, this NKDE is discontinuous, which is counter-intuitive, and leads to sharp differences between density values in the network, and @@ -299,8 +299,8 @@

The continuous method verbose = FALSE)
 d3_plot_situation(all_lines, event, samples_pts, continuous_kernel, scales = c(1,1,3))
-
-

As one can see, the values of the NKDE are continuous, and the +

+

As one can see, the values of the NKDE are continuous, and the density values close to the events have been adjusted. This method produces smoother results than the discontinuous method.

what to retain about the discontinuous method

@@ -337,8 +337,8 @@

The continuous method verbose = FALSE)
 d3_plot_situation(all_lines, events, samples_pts, continuous_kernel, scales = c(1,1,3))
-
-

The larger the bandwidth, the more common links the events share in +

+

The larger the bandwidth, the more common links the events share in their range.

 
@@ -350,8 +350,8 @@ 

The continuous method verbose = FALSE) d3_plot_situation(all_lines, events, samples_pts, continuous_kernel, scales = c(1,1,3))

-
- +
+

References diff --git a/docs/articles/web_vignettes/SpaceTimeDBscan.html b/docs/articles/web_vignettes/SpaceTimeDBscan.html index a90cb1aa..99ef5f49 100644 --- a/docs/articles/web_vignettes/SpaceTimeDBscan.html +++ b/docs/articles/web_vignettes/SpaceTimeDBscan.html @@ -110,7 +110,7 @@

Spatio-temporal dbscan with network

Jeremy Gelb

-

2023-11-01

+

2024-01-01

Source: vignettes/web_vignettes/SpaceTimeDBscan.Rmd @@ -213,8 +213,8 @@

Analysing the results tm_dots("black", alpha = 0.3, size = 0.01) + tm_shape(in_cluster) + tm_dots("cluster")

-
-

And finally, we will plot the time period of each cluster.

+
+

And finally, we will plot the time period of each cluster.

 
 ggplot(in_cluster) + 
diff --git a/docs/news/index.html b/docs/news/index.html
index 69aa6757..97f320d3 100644
--- a/docs/news/index.html
+++ b/docs/news/index.html
@@ -96,6 +96,7 @@ 

Changes in functions

corrected bugs

When using ESC-NKDE with bandwidths selection, negative density values could be obtained (with a small parameter of depth). It led to NAN values in LOO score. They are now treated like zeros.

+

A small bug has been corrected in the functions network_knn and network_knn.mc caused with grids with very few number of points.

diff --git a/docs/pkgdown.yml b/docs/pkgdown.yml index 316258ee..1adb8b91 100644 --- a/docs/pkgdown.yml +++ b/docs/pkgdown.yml @@ -11,5 +11,5 @@ articles: AdaptiveBW: web_vignettes/AdaptiveBW.html NKDEdetailed: web_vignettes/NKDEdetailed.html SpaceTimeDBscan: web_vignettes/SpaceTimeDBscan.html -last_built: 2023-11-02T00:46Z +last_built: 2024-01-02T03:23Z diff --git a/docs/reference/figures/unnamed-chunk-8-1.png b/docs/reference/figures/unnamed-chunk-8-1.png index 4e9ac0bf..219ffe2c 100644 Binary files a/docs/reference/figures/unnamed-chunk-8-1.png and b/docs/reference/figures/unnamed-chunk-8-1.png differ diff --git a/inst/extdata/results_vignette_kfunc.rda b/inst/extdata/results_vignette_kfunc.rda index 7d325a65..27f2dc1e 100644 Binary files a/inst/extdata/results_vignette_kfunc.rda and b/inst/extdata/results_vignette_kfunc.rda differ diff --git a/inst/extdata/results_vignette_network_build.rda b/inst/extdata/results_vignette_network_build.rda index bedf37b5..27356482 100644 Binary files a/inst/extdata/results_vignette_network_build.rda and b/inst/extdata/results_vignette_network_build.rda differ diff --git a/inst/extdata/results_vignette_nkde.rda b/inst/extdata/results_vignette_nkde.rda index d4796d9e..a5f55d6c 100644 Binary files a/inst/extdata/results_vignette_nkde.rda and b/inst/extdata/results_vignette_nkde.rda differ diff --git a/inst/extdata/results_vignette_tnkde.rda b/inst/extdata/results_vignette_tnkde.rda index f18647b5..a1dbc870 100644 Binary files a/inst/extdata/results_vignette_tnkde.rda and b/inst/extdata/results_vignette_tnkde.rda differ diff --git a/inst/extdata/results_vignette_wmat.rda b/inst/extdata/results_vignette_wmat.rda index 87b73aa9..adc514da 100644 Binary files a/inst/extdata/results_vignette_wmat.rda and b/inst/extdata/results_vignette_wmat.rda differ diff --git a/man/figures/unnamed-chunk-8-1.png b/man/figures/unnamed-chunk-8-1.png index 4e9ac0bf..219ffe2c 100644 Binary files a/man/figures/unnamed-chunk-8-1.png and b/man/figures/unnamed-chunk-8-1.png differ diff --git a/test.log b/test.log new file mode 100644 index 00000000..e77eef18 --- /dev/null +++ b/test.log @@ -0,0 +1,297 @@ +✔ | F W S OK | Context + ⠏ | 0 | adaptive_bw_sf ⠏ | 0 | testing functions for adptative bandwidth in NKDE [1] "checking inputs ..." +[1] "prior data preparation ..." +[1] "start calculating the kernel values for the adaptive bandwidth..." +[1] " quadra 1/3" +[1] " build graph ..." +[1] " calculating NKDE values ..." +[1] " quadra 2/3" +[1] " build graph ..." +[1] " calculating NKDE values ..." +[1] " quadra 3/3" +[1] " build graph ..." +[1] " calculating NKDE values ..." +[1] "Splitting the data with the spatial grid ..." +[1] "start calculating the kernel values ..." +[1] " quadra 1/3" +[1] " build graph ..." +[1] " calculating NKDE values ..." +[1] " quadra 2/3" +[1] " build graph ..." +[1] " calculating NKDE values ..." +[1] " quadra 3/3" +[1] " build graph ..." +[1] " calculating NKDE values ..." +[1] "combining the results ..." + ⠋ | 1 | testing functions for adptative bandwidth in NKDE [1] "precalculated bws : 3.73661539510889;2.93822642599985;2.45923541821083" +[1] "checking inputs ..." +[1] "prior data preparation ..." +[1] "calculating the local bandwidth ..." +[1] "start calculating the kernel values for the adaptive bandwidth..." +[1] "starting the kernel calculations ..." +[1] "starting the kernel calculations ..." +[1] "combining the results ..." + ⠙ | 2 | testing functions for adptative bandwidth in NKDE [1] "checking inputs ..." +[1] "prior data preparation ..." +[1] "start calculating the kernel values for the adaptive bandwidth..." +[1] " quadra 1/3" +[1] " build graph ..." +[1] " calculating NKDE values ..." +[1] " quadra 2/3" +[1] " build graph ..." +[1] " calculating NKDE values ..." +[1] " quadra 3/3" +[1] " build graph ..." +[1] " calculating NKDE values ..." +[1] "Splitting the data with the spatial grid ..." +[1] "start calculating the kernel values ..." +[1] " quadra 1/3" +[1] " build graph ..." +[1] " calculating NKDE values ..." +[1] " quadra 2/3" +[1] " build graph ..." +[1] " calculating NKDE values ..." +[1] " quadra 3/3" +[1] " build graph ..." +[1] " calculating NKDE values ..." +[1] "combining the results ..." + ⠹ | 3 | testing functions for adptative bandwidth in NKDE ✔ | 3 | testing functions for adptative bandwidth in NKDE [1.7s] + ⠏ | 0 | adaptive_bw_tnkde_interaction_sf ⠏ | 0 | testing functions for adptative bandwidth with interaction in TNKDE [1] "checking inputs ..." +[1] "prior data preparation ..." +[1] "start calculating the kernel values for the adaptive bandwidth..." +[1] "start calculating the kernel values ..." +[1] " quadra 1/1" +[1] " build graph ..." +[1] " calculating NKDE values ..." +[1] "combining the results ..." + ⠋ | 1 | testing functions for adptative bandwidth with interaction in TNKDE [1] "checking inputs ..." +[1] "prior data preparation ..." +[1] "start calculating the kernel values for the adaptive bandwidth..." +[1] "starting the kernel calculations ..." +[1] " build graph ..." +[1] " calculating NKDE values ..." +[1] "combining the results ..." + ⠙ | 2 | testing functions for adptative bandwidth with interaction in TNKDE [1] "start calculating the kernel values for the adaptive bandwidth..." + ⠹ | 3 | testing functions for adptative bandwidth with interaction in TNKDE [1] "start calculating the kernel values for the adaptive bandwidth..." + ⠦ | 7 | testing functions for adptative bandwidth with interaction in TNKDE ✔ | 8 | testing functions for adptative bandwidth with interaction in TNKDE [1.3s] + ⠏ | 0 | base_kernels_sf ⠏ | 0 | testing the base kernel implemented ✔ | 2 | testing the base kernel implemented + ⠏ | 0 | basics_listw_sf ⠏ | 0 | basic tests for the neighbour matrices creation [1] "generating the starting points" +[1] "snapping the points to the lines (only once)" +[1] "starting the network part" +[1] "working on quadra : 1/1" +[1] "adding the points as vertices to nearest lines" +[1] "splitting the lines for the network" +[1] "generating the network" +[1] "calculating the distances on the network" +[1] "building the final listw object" +[1] "finally generating the listw object ..." + ⠙ | 2 | basic tests for the neighbour matrices creation ✔ | 2 | basic tests for the neighbour matrices creation [0.3s] + ⠏ | 0 | border_correction_sf ⠏ | 0 | testing function used to apply border correction in NKDE ⠙ | 2 | testing function used to apply border correction in NKDE ⠹ | 3 | testing function used to apply border correction in NKDE ⠸ | 4 | testing function used to apply border correction in NKDE ✔ | 4 | testing function used to apply border correction in NKDE [1.1s] + ⠏ | 0 | bw_selection_sf ⠏ | 0 | testing the bandwidth selection functions [1] 1 +[1] 1 + ⠙ | 2 | testing the bandwidth selection functions [1] 1 +[1] 1 + ⠸ | 4 | testing the bandwidth selection functions [1] "checking inputs ..." +[1] "prior data preparation ..." +[1] "Calculating the correction factor if required" +[1] "Splittig the dataset within the grid ..." +[1] "start calculating the CV values ..." + | | | 0%[1] 1 + | |=====================================================================================================================| 100% ⠴ | 6 | testing the bandwidth selection functions [1] "checking inputs ..." +[1] "prior data preparation ..." +[1] "Calculating the correction factor if required" +[1] "start calculating the CV values ..." + | | | 0% ⠧ | 8 | testing the bandwidth selection functions ⠇ | 9 | testing the bandwidth selection functions [1] "checking inputs ..." +[1] "prior data preparation ..." +[1] "Calculating the correction factor if required" +[1] "start calculating the CV values ..." +[1] "checking inputs ..." +[1] "prior data preparation ..." +[1] "Calculating the correction factor if required" +[1] "Splittig the dataset within the grid ..." +[1] "start calculating the CV values ..." + | | | 0%[1] 1 + | |============================= | 25%[1] 2 + | |========================================================== | 50%[1] 3 + | |======================================================================================== | 75%[1] 4 + | |=====================================================================================================================| 100% ⠏ | 10 | testing the bandwidth selection functions [1] "checking inputs ..." +[1] "prior data preparation ..." +[1] "Calculating the correction factor if required" +[1] "start calculating the CV values ..." + bw cvl_scores +1 200 3.068423e+12 +2 300 9.065686e+12 +3 400 1.812569e+13 + bw cvl_scores +1 200 3.068423e+12 +2 300 9.065686e+12 +3 400 1.812569e+13 + bw cvl_scores +1 200 3.068423e+12 +2 300 9.065686e+12 +3 400 1.812569e+13 + ⠋ | 11 | testing the bandwidth selection functions [1] "checking inputs ..." +[1] "prior data preparation ..." +[1] "Calculating the correction factor if required" +[1] "start calculating the kernel values for the adaptive bandwidth..." +[1] " quadra 1/3" +[1] " build graph ..." +[1] " calculating NKDE values ..." +[1] " quadra 2/3" +[1] " build graph ..." +[1] " calculating NKDE values ..." +[1] " quadra 3/3" +[1] " build graph ..." +[1] " calculating NKDE values ..." +[1] "Splittig the dataset within the grid ..." +[1] "start calculating the CV values ..." + | | | 0%[1] 1 + | |======================================= | 33%[1] 2 + | |============================================================================== | 67%[1] 3 + | |=====================================================================================================================| 100% ⠙ | 12 | testing the bandwidth selection functions [1] 1 + ⠹ | 13 | testing the bandwidth selection functions ✔ | 13 | testing the bandwidth selection functions [7.6s] + ⠏ | 0 | bw_selection_tnkde_sf ⠏ | 0 | testing the bandwidth selection functions for tnkde ⠋ | 1 | testing the bandwidth selection functions for tnkde ⠙ | 2 | testing the bandwidth selection functions for tnkde [1] "prior data preparation ..." +[1] "Calculating the correction factor if required" +[1] "start calculating the CV values ..." + | | | 0% | |=====================================================================================================================| 100%[1] "prior data preparation ..." +[1] "Calculating the correction factor if required" +[1] "start calculating the CV values ..." + | | | 0% | |=====================================================================================================================| 100% ⠸ | 4 | testing the bandwidth selection functions for tnkde [1] "prior data preparation ..." +[1] "Calculating the correction factor if required" +[1] "start calculating the CV values ..." + | | | 0% | |=====================================================================================================================| 100%[1] "prior data preparation ..." +[1] "Calculating the correction factor if required" +[1] "splitting the data by the grid..." +[1] "start calculating the CV values ..." + ⠼ | 5 | testing the bandwidth selection functions for tnkde ⠴ | 6 | testing the bandwidth selection functions for tnkde ⠦ | 7 | testing the bandwidth selection functions for tnkde ⠧ | 8 | testing the bandwidth selection functions for tnkde ✔ | 8 | testing the bandwidth selection functions for tnkde [1.5s] + ⠏ | 0 | geometrical_functions_sf ⠏ | 0 | Testing the geometrical functions of the package | | | 0% | |======================= | 20% | |=============================================== | 40% | |====================================================================== | 60% | |============================================================================================== | 80% | |=====================================================================================================================| 100% ⠴ | 6 | Testing the geometrical functions of the package ⠋ | 11 | Testing the geometrical functions of the package | | | 0% | |============================= | 25% | |========================================================== | 50% | |======================================================================================== | 75% | |=====================================================================================================================| 100% ✔ | 12 | Testing the geometrical functions of the package [0.3s] + ⠏ | 0 | graph_functions_sf ⠏ | 0 | testing the graph functions ✔ | 2 | testing the graph functions + ⠏ | 0 | isochrone_sf ⠏ | 0 | testing the isochrone creation functions ⠙ | 2 | testing the isochrone creation functions ✔ | 2 | testing the isochrone creation functions [0.2s] + ⠏ | 0 | testing the isochrone creation functions ⠋ | 1 | testing the isochrone creation functions ✔ | 1 | testing the isochrone creation functions [0.1s] + ⠏ | 0 | k_functions_sf ⠏ | 0 | testing functions for k and g function analysis [1] "Preparing data ..." +[1] "Snapping points on lines ..." +[1] "Building graph ..." +[1] "Calculating k and g functions ..." +[1] "Calculating the simulations ..." + | | | 0% | |======================= | 20% | |=============================================== | 40% | |====================================================================== | 60% | |============================================================================================== | 80% | |=====================================================================================================================| 100% ⠋ | 1 0 | testing functions for k and g function analysis [1] "Preparing data ..." +[1] "Snapping points on lines ..." +[1] "Building graph ..." +[1] "Calculating k and g functions ..." +[1] "Calculating the simulations ..." + | | | 0% | |======================= | 20% | |=============================================== | 40% | |====================================================================== | 60% | |============================================================================================== | 80% | |=====================================================================================================================| 100% ⠹ | 1 2 | testing functions for k and g function analysis ⠼ | 1 4 | testing functions for k and g function analysis [1] "Preparing data ..." +[1] "Snapping points on lines ..." +[1] "Building graph ..." +[1] "Calculating k and g functions ..." +[1] "Calculating the simulations ..." + ⠴ | 2 4 | testing functions for k and g function analysis [1] "Preparing data ..." +[1] "Snapping points on lines ..." +[1] "Building graph ..." +[1] "Calculating k and g functions ..." +[1] "Calculating the simulations ..." + ⠧ | 3 5 | testing functions for k and g function analysis [1] "Preparing data ..." +[1] "Snapping points on lines ..." +[1] "Building graph ..." +[1] "Calculating k and g functions ..." +[1] "Calculating the simulations ..." + | | | 0% ⠏ | 4 6 | testing functions for k and g function analysis [1] "Preparing data ..." +[1] "Snapping points on lines ..." +[1] "Building graph ..." +[1] "Calculating k and g functions ..." +[1] "Calculating the simulations ..." + | | | 0% | |======================= | 20% | |=============================================== | 40% | |====================================================================== | 60% | |============================================================================================== | 80% | |=====================================================================================================================| 100% ⠙ | 4 8 | testing functions for k and g function analysis ✔ | 4 8 | testing functions for k and g function analysis [1.4s] +─────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────── +Warning ('test_k_functions_sf.R:28:3'): Testing the simple k function +`aes_string()` was deprecated in ggplot2 3.0.0. +ℹ Please use tidy evaluation idioms with `aes()`. +ℹ See also `vignette("ggplot2-in-packages")` for more information. +Backtrace: + 1. spNetwork::kfunctions(...) + at test_k_functions_sf.R:28:2 + 4. ggplot2::aes_string(x = "distances", ymin = "lower_k", ymax = "upper_k") + 5. ggplot2:::deprecate_soft0(...) + +Warning ('test_k_functions_sf.R:332:3'): Testing the multicore simple k function +UNRELIABLE VALUE: One of the ‘future.apply’ iterations (‘future_lapply-1’) unexpectedly generated random numbers without declaring so. There is a risk that those random numbers are not statistically sound and the overall results might be invalid. To fix this, specify 'future.seed=TRUE'. This ensures that proper, parallel-safe random numbers are produced via the L'Ecuyer-CMRG method. To disable this check, use 'future.seed = NULL', or set option 'future.rng.onMisuse' to "ignore". +Backtrace: + 1. spNetwork::kfunctions.mc(...) + at test_k_functions_sf.R:332:2 + 2. progressr::with_progress(...) + at spNetwork/R/k_functions_sf.R:599:6 + 3. progressr:::flush_conditions(conditions) + at spNetwork/R/k_functions_sf.R:599:6 + +Warning ('test_k_functions_sf.R:368:3'): Testing the multicore cross k function +UNRELIABLE VALUE: One of the ‘future.apply’ iterations (‘future_lapply-1’) unexpectedly generated random numbers without declaring so. There is a risk that those random numbers are not statistically sound and the overall results might be invalid. To fix this, specify 'future.seed=TRUE'. This ensures that proper, parallel-safe random numbers are produced via the L'Ecuyer-CMRG method. To disable this check, use 'future.seed = NULL', or set option 'future.rng.onMisuse' to "ignore". +Backtrace: + 1. spNetwork::cross_kfunctions.mc(...) + at test_k_functions_sf.R:368:2 + 2. progressr::with_progress(...) + at spNetwork/R/k_functions_sf.R:1077:4 + 3. progressr:::flush_conditions(conditions) + at spNetwork/R/k_functions_sf.R:1077:4 + +Warning ('test_k_functions_sf.R:407:3'): Testing the simple k function in space-time +UNRELIABLE VALUE: One of the ‘future.apply’ iterations (‘future_lapply-1’) unexpectedly generated random numbers without declaring so. There is a risk that those random numbers are not statistically sound and the overall results might be invalid. To fix this, specify 'future.seed=TRUE'. This ensures that proper, parallel-safe random numbers are produced via the L'Ecuyer-CMRG method. To disable this check, use 'future.seed = NULL', or set option 'future.rng.onMisuse' to "ignore". +Backtrace: + 1. spNetwork::k_nt_functions.mc(...) + at test_k_functions_sf.R:407:2 + 2. progressr::with_progress(...) + at spNetwork/R/k_functions_nettime_sf.R:452:6 + 3. progressr:::flush_conditions(conditions) + at spNetwork/R/k_functions_nettime_sf.R:452:6 +─────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────── + ⠏ | 0 | kernels ⠏ | 0 | testing the kernel functions [1] "checking inputs ..." +[1] "prior data preparation ..." +[1] "Splitting the data with the spatial grid ..." +[1] "start calculating the kernel values ..." +[1] " quadra 1/1" +[1] " build graph ..." +[1] " calculating NKDE values ..." + | | | 0% | |=====================================================================================================================| 100%[1] "combining the results ..." + ⠙ | 2 | testing the kernel functions ⠸ | 4 | testing the kernel functions ⠼ | 5 | testing the kernel functions ⠴ | 6 | testing the kernel functions ✔ | 6 | testing the kernel functions [0.6s] + ⠏ | 0 | knn_sf ⠏ | 0 | testing the knn functions ⠙ | 2 | testing the knn functions ⠸ | 4 | testing the knn functions ✔ | 5 | testing the knn functions [0.4s] + ⠏ | 0 | listw_symetry_sf ⠏ | 0 | symetry and comparison between simple and multicore ⠋ | 1 | symetry and comparison between simple and multicore ✔ | 1 | symetry and comparison between simple and multicore [7.3s] + ⠏ | 0 | nkde_execution_functions_sf ⠏ | 0 | testing function used to perform NKDE ✔ | 1 | testing function used to perform NKDE + ⠏ | 0 | nkde_manual_bws ⠏ | 0 | testing the kernel functions with manually changed bandwidths ⠙ | 2 | testing the kernel functions with manually changed bandwidths ✔ | 2 | testing the kernel functions with manually changed bandwidths [0.2s] + ⠏ | 0 | spatial_indexing ⠏ | 0 | testing the spatial indexing functions ✔ | 1 | testing the spatial indexing functions + ⠏ | 0 | splitting_functions ⠏ | 0 | testing the base kernel implemented ⠋ | 1 | testing the base kernel implemented ⠙ | 2 | testing the base kernel implemented ✔ | 2 | testing the base kernel implemented [0.8s] + ⠏ | 0 | tnkde_sf ⠏ | 0 | testing the kernel functions [1] "checking inputs ..." +[1] "prior data preparation ..." +[1] "start calculating the kernel values ..." +[1] " quadra 1/1" +[1] " build graph ..." +[1] " calculating NKDE values ..." + | | | 0% | |=====================================================================================================================| 100%[1] "combining the results ..." +[1] "checking inputs ..." +[1] "prior data preparation ..." +[1] "start calculating the kernel values ..." +[1] " quadra 1/1" +[1] " build graph ..." +[1] " calculating NKDE values ..." +[1] "combining the results ..." + ⠙ | 2 | testing the kernel functions [1] "checking inputs ..." +[1] "prior data preparation ..." +[1] "start calculating the kernel values ..." +[1] " quadra 1/1" +[1] " build graph ..." +[1] " calculating NKDE values ..." +[1] "combining the results ..." + ⠹ | 3 | testing the kernel functions [1] "checking inputs ..." +[1] "prior data preparation ..." +[1] "start calculating the kernel values ..." +[1] " quadra 1/1" +[1] " build graph ..." +[1] " calculating NKDE values ..." +[1] "combining the results ..." +[1] "checking inputs ..." +[1] "prior data preparation ..." +[1] "starting the kernel calculations ..." +[1] " build graph ..." +[1] " calculating NKDE values ..." +[1] "combining the results ..." + ⠸ | 4 | testing the kernel functions ✔ | 4 | testing the kernel functions [0.7s] + +══ Results ════════════════════════════════════════════════════════════════════════════════════════════════════════════════════ +Duration: 25.9 s + +[ FAIL 0 | WARN 4 | SKIP 0 | PASS 87 ] diff --git a/vignettes/NKDE.Rmd b/vignettes/NKDE.Rmd index b167f642..df46d813 100644 --- a/vignettes/NKDE.Rmd +++ b/vignettes/NKDE.Rmd @@ -576,24 +576,24 @@ For both, one can use the multicore version of the functions to speed up the cal ```{r eval = FALSE} bws_selection_cv <- bw_cv_likelihood_calc( - bws = seq(200,700,50), + bws = seq(50,700,50), lines = mtl_network, events = bike_accidents, w = rep(1,nrow(bike_accidents)), kernel_name = "quartic", method = "discontinuous", diggle_correction = FALSE, study_area = NULL, max_depth = 8, - digits=2, tol=0.1, agg=5, + digits=2, tol=0.1, agg=5, sparse=TRUE, grid_shape=c(1,1), verbose=FALSE, check=TRUE) bws_selection_cvl <- bw_cvl_calc( - bws = seq(50,700, 50), + bws = seq(50,700,50), lines = mtl_network, events = bike_accidents, w = rep(1,nrow(bike_accidents)), kernel_name = "quartic", method = "discontinuous", diggle_correction = FALSE, study_area = NULL, max_depth = 8, - digits=2, tol=0.1, agg=5, + digits=2, tol=0.1, agg=5, sparse=TRUE, grid_shape=c(1,1), verbose=FALSE, check=TRUE) @@ -602,22 +602,52 @@ bws_selection_cvl <- bw_cvl_calc( It is also possible to do bandwidth selection with adaptive bandwidth. In that context, we search for the optimal global bandwidth. It requires more computation because for each bandwidth, we need to estimate the density twice, one with the global bandwidth, and one with the local bandwidths obtained from the density estimation. ```{r eval = FALSE} -bws_selection_cv_adpt <- bw_cv_likelihood_calc( - bws = seq(200,700,50), - lines = mtl_network, events = bike_accidents, +bws_selection_cv_adpt_dis <- bw_cv_likelihood_calc( + bws = seq(50,700,50), + trim_bws = seq(50,700,50)*2, + lines = mtl_network, + events = bike_accidents, w = rep(1,nrow(bike_accidents)), adaptive = TRUE, - kernel_name = "quartic", method = "discontinuous", - diggle_correction = FALSE, study_area = NULL, + kernel_name = "quartic", + method = "discontinuous", + diggle_correction = FALSE, + study_area = NULL, max_depth = 8, - digits=2, tol=0.1, agg=5, - sparse=TRUE, grid_shape=c(1,1), - verbose=FALSE, check=TRUE) + digits=2, + tol=0.1, + agg=5, + sparse=TRUE, + grid_shape=c(1,1), + verbose=TRUE, + check=TRUE) + +bws_selection_cv_adpt_cont <- bw_cv_likelihood_calc( + bws = seq(50,700,50), + trim_bws = seq(50,700,50)*2, + lines = mtl_network, + events = bike_accidents, + w = rep(1,nrow(bike_accidents)), + adaptive = TRUE, + kernel_name = "quartic", + method = "continuous", + diggle_correction = FALSE, + study_area = NULL, + max_depth = 8, + digits=2, + tol=0.1, + agg=5, + sparse=TRUE, + grid_shape=c(1,1), + verbose=TRUE, + check=TRUE) cv_values <- data.frame( "bw" = bws_selection_cv$bw, "cv_likelihood" = bws_selection_cv$cv_scores, - "cvl_crit" = bws_selection_cvl$cvl_scores + "cvl_crit" = bws_selection_cvl$cvl_scores, + "cv_likelihood_adpt_dis" = bws_selection_cv_adpt_dis$cv_scores, + "cv_likelihood_adpt_cont" = bws_selection_cv_adpt_cont$cv_scores ) ``` @@ -625,7 +655,7 @@ cv_values <- data.frame( knitr::kable(cv_values, digits = 2,row.names = FALSE) ``` -As one can see here, the cross-validation approach suggests a value between 550 and 600 metres, while the Cronie and Van Lieshout's Criterion would suggest a bandwidth inferior to 50 metres and is probably not adapted for this dataset. With an adaptive bandwidth, the best bandwidth would be around 500 and 550 metres. +As one can see here, the cross-validation approach suggests a value between 550 and 600 metres, while the Cronie and Van Lieshout's Criterion would suggest a bandwidth inferior to 50 metres and is probably not adapted for this dataset. With an adaptive bandwidth, the best bandwidth would be around 450 and 500 metres for the discontinuous kernel, and between 550 and 600 for the continuous kernel. ## Complementary functions diff --git a/vignettes/TNKDE.Rmd b/vignettes/TNKDE.Rmd index a7332249..835c2a52 100644 --- a/vignettes/TNKDE.Rmd +++ b/vignettes/TNKDE.Rmd @@ -291,17 +291,15 @@ One can decide to use adaptive bandwidth for the TNKDE. They are calculated usin We can also select the global bandwidth by cross validation with the function `bw_tnkde_cv_likelihood_calc`. ```{r message=FALSE, warning=FALSE, eval = FALSE} -cv_scores_adpt <- bw_tnkde_cv_likelihood_calc( - bw_net_range = c(200,1100), - bw_net_step = 100, - bw_time_range = c(10,70), - bw_time_step = 10, +cv_scores_adpt_dis <- bw_tnkde_cv_likelihood_calc( + bws_net = seq(200,1100,100), + bws_time = seq(10,70,10), lines = mtl_network, events = bike_accidents, time_field = "Time", adaptive = TRUE, - trim_net_bws = seq(200,1100,100), - trim_time_bws = seq(10,70,10), + trim_net_bws = seq(200,1100,100)*2, + trim_time_bws = seq(10,70,10)*2, w = rep(1, nrow(bike_accidents)), kernel_name = "quartic", method = "discontinuous", @@ -316,12 +314,36 @@ cv_scores_adpt <- bw_tnkde_cv_likelihood_calc( sub_sample=1, verbose = FALSE, check = TRUE) + +cv_scores_adpt_cont <- bw_tnkde_cv_likelihood_calc.mc( + bws_net = seq(200,1100,100), + bws_time = seq(10,70,10), + lines = mtl_network, + events = bike_accidents, + time_field = "Time", + adaptive = TRUE, + trim_net_bws = seq(200,1100,100)*2, + trim_time_bws = seq(10,70,10)*2, + w = rep(1, nrow(bike_accidents)), + kernel_name = "quartic", + method = "continuous", + diggle_correction = FALSE, + study_area = NULL, + max_depth = 10, + digits = 2, + tol = 0.1, + agg = 15, + sparse=TRUE, + grid_shape=c(2,2), + sub_sample=1, + verbose = FALSE, + check = TRUE) ``` ```{r message=FALSE, warning=FALSE} -knitr::kable(cv_scores_adpt) +knitr::kable(cv_scores_adpt_cont) ``` -The best score is reached with a network and temporal bandwidths of 800m and 60 days. +For a discontinuous NKDE The best score is reached with a network and temporal bandwidths of 800m and 60 days. ```{r message=FALSE, warning=FALSE, eval = FALSE} tnkde_densities <- tnkde(lines = mtl_network, diff --git a/vignettes/images/animated_map.gif b/vignettes/images/animated_map.gif index a85183c1..81825931 100644 Binary files a/vignettes/images/animated_map.gif and b/vignettes/images/animated_map.gif differ diff --git a/vignettes/web_vignettes/AdaptiveBW.Rmd b/vignettes/web_vignettes/AdaptiveBW.Rmd index 520a7867..b2dda56d 100644 --- a/vignettes/web_vignettes/AdaptiveBW.Rmd +++ b/vignettes/web_vignettes/AdaptiveBW.Rmd @@ -98,6 +98,7 @@ densities <- nkde(lines = mtl_network, ``` We can now map the result + ```{r message=FALSE, warning=FALSE, fig.align= "center", fig.width = 8, fig.height = 8} # rescaling to help the mapping @@ -181,6 +182,7 @@ tmap_animation(all_maps, filename = "images/animated_densities.gif", ```{r message=FALSE, warning=FALSE, fig.align= "center", out.width = "100%"} knitr::include_graphics("images/animated_densities.gif") ``` + With no suprise, the larger the bandwidth, the smoother the results ## Finding the optimal bandwidth @@ -189,12 +191,13 @@ To complete the visual inspection we did above, we could select an optimal set o To do so, we can use the function `bw_cv_likelihood_calc` and we must create a matrix of bandwidths that will be evaluated. We will compare below the scores obtained for the local bandwidths based on the 5th to 25th nearest neighbours. + ```{r message=FALSE, warning=FALSE} nearest_nb_mat <- knn_dists$distances[,5:25] colnames(nearest_nb_mat) <- 5:25 -bw_scores <- bw_cv_likelihood_calc(lines = mtl_network, +bw_scores_dis <- bw_cv_likelihood_calc(lines = mtl_network, events = bike_accidents_agg, w = bike_accidents_agg$weight, kernel_name = 'quartic', @@ -202,21 +205,41 @@ bw_scores <- bw_cv_likelihood_calc(lines = mtl_network, adaptive = TRUE, mat_bws = nearest_nb_mat, digits = 1, - tol = 0.01, - max_depth = 10, agg = NULL, sparse = TRUE) + tol = 0.01, + grid_shape = c(2,2), + max_depth = 8, agg = NULL, sparse = TRUE) + +bw_scores_cont <- bw_cv_likelihood_calc(lines = mtl_network, + events = bike_accidents_agg, + w = bike_accidents_agg$weight, + kernel_name = 'quartic', + method = 'continuous', + adaptive = TRUE, + mat_bws = nearest_nb_mat, + digits = 1, + tol = 0.01, + grid_shape = c(2,2), + max_depth = 8, agg = NULL, sparse = TRUE) df <- data.frame( k = 5:25, - score = bw_scores[,2] + score_dis = bw_scores_dis[,2], + score_cont = bw_scores_cont[,2] ) + ggplot(df) + - geom_line(aes(x = k, y = score)) + - geom_point(aes(x = k, y = score), color = 'red') + + geom_line(aes(x = k, y = score_dis)) + + geom_point(aes(x = k, y = score_dis, color = 'discontinuous')) + + geom_line(aes(x = k, y = score_cont)) + + geom_point(aes(x = k, y = score_cont, color = 'continuous')) + + scale_color_manual(values = c("discontinuous" = "red", + "continuous" = "blue")) + + labs(color = '') + theme_bw() ``` -Based on these results, we decide to use the 8th nearest neighbour as an optimal bandwidth. It clearly minimizes the likelihood cross validation criterion. +Based on these results, we decide to use the 8th nearest neighbour as an optimal bandwidth. It clearly maximizes the likelihood cross validation criterion for both the continuous and the discontinuous NKDE. # K nearest neighbours adaptive bandwidth for TNKDE @@ -364,7 +387,7 @@ for(j in 5:25){ } -values <- bw_tnkde_cv_likelihood_calc(lines = mtl_network, +scores_dis <- bw_tnkde_cv_likelihood_calc(lines = mtl_network, events = bike_accidents, time_field = "int_time", w = rep(1,nrow(bike_accidents)), @@ -388,16 +411,48 @@ values <- bw_tnkde_cv_likelihood_calc(lines = mtl_network, verbose=FALSE, check=TRUE) +scores_cont <- bw_tnkde_cv_likelihood_calc(lines = mtl_network, + events = bike_accidents, + time_field = "int_time", + w = rep(1,nrow(bike_accidents)), + kernel_name = 'quartic', + method = 'continuous', + arr_bws_net = network_bws, + arr_bws_time = time_bws, + diggle_correction = FALSE, + study_area = NULL, + adaptive = TRUE, + trim_net_bws = NULL, + trim_time_bws = NULL, + max_depth = 10, + digits=5, + tol=0.1, + agg=NULL, + sparse=TRUE, + zero_strat = "min_double", + grid_shape=c(1,1), + sub_sample=1, + verbose=FALSE, + check=TRUE) + df <- data.frame( k = 5:25, - score = values[,1] + score_dis = scores_dis[,1], + score_cont = scores_cont[,1] ) + ggplot(df) + - geom_line(aes(x = k, y = score)) + - geom_point(aes(x = k, y = score), color = 'red') + + geom_line(aes(x = k, y = score_dis)) + + geom_point(aes(x = k, y = score_dis, color = 'discontinuous')) + + geom_line(aes(x = k, y = score_cont)) + + geom_point(aes(x = k, y = score_cont, color = 'continuous')) + + scale_color_manual(values = c("discontinuous" = "red", + "continuous" = "blue")) + + labs(color = '') + theme_bw() + ``` -By checking the table above, we decide to select the 10th nearest neighbour because it clearly corresponds to the major elbow in the above figure. Also, the 21th neighbour gives the highest score. +By checking the table above, we could decide to select the 10th nearest neighbour because it clearly corresponds to the major elbow in the above figure. Also, the 21th neighbour gives the highest score.