Merge branch 'develop' into check-on-newer-R

docker/depends/pecan_package_dependencies.csv

@@ -474,7 +474,6 @@
 "RPostgres","*","base/db","Suggests",FALSE
 "RPostgreSQL","*","base/db","Suggests",FALSE
 "RPostgreSQL","*","models/biocro","Suggests",FALSE
-"RPostgreSQL","*","modules/data.land","Suggests",FALSE
 "Rpreles","*","models/preles","Suggests",FALSE
 "RSQLite","*","base/db","Suggests",FALSE
 "sf","*","modules/assim.sequential","Suggests",FALSE

modules/data.land/DESCRIPTION

@@ -61,7 +61,6 @@ Suggests:
     PEcAn.settings,
     redland,
     raster,
-    RPostgreSQL,
     testthat (>= 1.0.2)
 License: BSD_3_clause + file LICENSE
 Copyright: Authors

modules/data.land/R/IC_BADM_Utilities.R

@@ -21,7 +21,7 @@ Read.IC.info.BADM <-function(lat, long){
   #Reading in the DB
   #
   U.S.SB <-
-    read.csv(system.file("data","BADM.csv", package = "PEcAn.data.land"))
+    utils::read.csv(system.file("data","BADM.csv", package = "PEcAn.data.land"))
 
 
   Regions <- EPA_ecoregion_finder(lat, long)
@@ -31,8 +31,8 @@ Read.IC.info.BADM <-function(lat, long){
   #L2
   biomass.df <- U.S.SB %>%
     dplyr::filter(
-      NA_L2CODE == Code_Level,
-      VARIABLE %>% grepl("ROOT_|AG_BIOMASS|SOIL_STOCK|SOIL_CHEM", .)
+      .data$NA_L2CODE == Code_Level,
+      grepl("ROOT_|AG_BIOMASS|SOIL_STOCK|SOIL_CHEM", .data$VARIABLE)
     ) %>%
     dplyr::select("SITE_ID", "GROUP_ID", "VARIABLE_GROUP", "VARIABLE", "DATAVALUE")
 
@@ -43,8 +43,8 @@ Read.IC.info.BADM <-function(lat, long){
 
     biomass.df <- U.S.SB %>%
       dplyr::filter(
-        NA_L1CODE == Code_Level,
-        VARIABLE %>% grepl("ROOT_|AG_BIOMASS|SOIL_STOCK|SOIL_CHEM", .)
+        .data$NA_L1CODE == Code_Level,
+        grepl("ROOT_|AG_BIOMASS|SOIL_STOCK|SOIL_CHEM", .data$VARIABLE)
       ) %>%
       dplyr::select("SITE_ID", "GROUP_ID", "VARIABLE_GROUP", "VARIABLE", "DATAVALUE")
   }
@@ -54,7 +54,7 @@ Read.IC.info.BADM <-function(lat, long){
   if (nrow(biomass.df) < 3)  {
     Code_Level <- "ALL"
     biomass.df <- U.S.SB %>%
-      dplyr::filter(VARIABLE %>% grepl("ROOT_|AG_BIOMASS|SOIL_STOCK|SOIL_CHEM", .)) %>%
+      dplyr::filter(grepl("ROOT_|AG_BIOMASS|SOIL_STOCK|SOIL_CHEM", .data$VARIABLE)) %>%
       dplyr::select("SITE_ID", "GROUP_ID", "VARIABLE_GROUP", "VARIABLE", "DATAVALUE")
   }
 
@@ -94,8 +94,8 @@ Read.IC.info.BADM <-function(lat, long){
           type <- type[-which(type == "*_BIOMASS")]
         #----------------- Unit conversion
         unit.in <- Gdf %>%
-          dplyr::filter(VARIABLE %>% grepl("UNIT", .)) %>%
-          dplyr::pull(DATAVALUE) %>% 
+          dplyr::filter(grepl("UNIT", .data$VARIABLE)) %>%
+          dplyr::pull(.data$DATAVALUE) %>% 
           as.character()
 
 
@@ -113,8 +113,8 @@ Read.IC.info.BADM <-function(lat, long){
           unit.ready <- "kg/m^2"
 
         Date.in <- Gdf %>%
-          dplyr::filter(VARIABLE %>% grepl("DATE", .)) %>%
-          dplyr::pull(DATAVALUE) %>%
+          dplyr::filter(grepl("DATE", .data$VARIABLE)) %>%
+          dplyr::pull(.data$DATAVALUE) %>%
           as.Date()
 
         if (length(Date.in) == 0)
@@ -124,8 +124,8 @@ Read.IC.info.BADM <-function(lat, long){
         # if it's biomass
         if (type == "*_BIOMASS") {
           Oregan.in <- Gdf %>%
-            dplyr::filter(VARIABLE %>% grepl("ORGAN", .)) %>%
-            dplyr::pull(DATAVALUE) 
+            dplyr::filter(grepl("ORGAN", .data$VARIABLE)) %>%
+            dplyr::pull(.data$DATAVALUE) 
 
 
           PlantWoodIni <-
@@ -134,8 +134,8 @@ Read.IC.info.BADM <-function(lat, long){
 
         } else if (type == "*SOIL") {
           val <- Gdf %>%
-            dplyr::filter(VARIABLE %>% grepl("SOIL_STOCK_C_ORG", .)) %>% #"SOIL_STOCK_C_ORG"
-            dplyr::pull(DATAVALUE) %>%
+            dplyr::filter(grepl("SOIL_STOCK_C_ORG", .data$VARIABLE)) %>%
+            dplyr::pull(.data$DATAVALUE) %>%
             as.numeric()
 
           if (length(val) > 0)

modules/data.land/R/InventoryGrowthFusion.R

@@ -404,8 +404,8 @@ model{
 
   ## JAGS initial conditions
   init   <- list()
-  if(is.mcmc.list(restart)){
-    init <- mcmc.list2init(restart)
+  if(coda::is.mcmc.list(restart)){
+    init <- PEcAn.utils::mcmc.list2init(restart)
     nchain <- length(init)
   } else {
     nchain <- 3
@@ -438,7 +438,7 @@ model{
   }
 
   PEcAn.logger::logger.info("RUN MCMC")
-  load.module("dic")
+  rjags::load.module("dic")
   for(k in avail.chunks){
 
     ## determine whether to sample states
@@ -457,17 +457,17 @@ model{
     save(jags.out,file=ofile)
 
     ## update restart
-    if(!is.null(restart) & ((is.logical(restart) && restart) || is.mcmc.list(restart))){
+    if(!is.null(restart) & ((is.logical(restart) && restart) || coda::is.mcmc.list(restart))){
       ofile <- paste("IGF",model,"RESTART.RData",sep=".")
-      jags.final <- coda.samples(model = j.model, variable.names = c("x",out.variables), n.iter = 1)
+      jags.final <- rjags::coda.samples(model = j.model, variable.names = c("x",out.variables), n.iter = 1)
       k_restart = k + 1  ## finished k, so would restart at k+1
       save(jags.final,k_restart,file=ofile)
     }
 
     ## check for convergence and break from loop early
-    D <- as.mcmc.list(lapply(jags.out,function(x){x[,'deviance']}))
+    D <- coda::as.mcmc.list(lapply(jags.out,function(x){x[,'deviance']}))
     gbr <- coda::gelman.diag(D)$psrf[1,1]
-    trend <- mean(sapply(D,function(x){coef(lm(x~seq_len(n.chunk)))[2]}))
+    trend <- mean(sapply(D,function(x){stats::coef(stats::lm(x~seq_len(n.chunk)))[2]}))
     if(gbr < 1.005 & abs(trend) < 0.5) break
   }
 

modules/data.land/R/InventoryGrowthFusionDiagnostics.R

@@ -9,13 +9,13 @@ InventoryGrowthFusionDiagnostics <- function(jags.out, combined=NULL) {
   out      <- as.matrix(jags.out)
   x.cols   <- which(substr(colnames(out), 1, 1) == "x")
   if(length(x.cols) > 0){
-    ci       <- apply(out[, x.cols], 2, quantile, c(0.025, 0.5, 0.975))
+    ci       <- apply(out[, x.cols], 2, stats::quantile, c(0.025, 0.5, 0.975))
     ci.names <- parse.MatrixNames(colnames(ci), numeric = TRUE)
 
     ### DBH par(mfrow=c(3,2))
     if(length(x.cols) > 0){
-      layout(matrix(1:8, 4, 2, byrow = TRUE))
-      ci       <- apply(out[, x.cols], 2, quantile, c(0.025, 0.5, 0.975))
+      graphics::layout(matrix(1:8, 4, 2, byrow = TRUE))
+      ci       <- apply(out[, x.cols], 2, stats::quantile, c(0.025, 0.5, 0.975))
       ci.names <- parse.MatrixNames(colnames(ci), numeric = TRUE)
 
       smp <- sample.int(data$ni, min(8, data$ni))
@@ -26,29 +26,29 @@ InventoryGrowthFusionDiagnostics <- function(jags.out, combined=NULL) {
         plot(data$time, ci[2, sel], type = "n", 
              ylim = range(rng), ylab = "DBH (cm)", main = i)
         PEcAn.visualization::ciEnvelope(data$time, ci[1, sel], ci[3, sel], col = "lightBlue")
-        points(data$time, data$z[i, ], pch = "+", cex = 1.5)
+        graphics::points(data$time, data$z[i, ], pch = "+", cex = 1.5)
         # lines(data$time,z0[i,],lty=2)
 
         ## growth
         sel      <- which(ci.names$row == i)
         inc.mcmc <- apply(out[, x.cols[sel]], 1, diff)
-        inc.ci   <- apply(inc.mcmc, 1, quantile, c(0.025, 0.5, 0.975)) * 5
+        inc.ci   <- apply(inc.mcmc, 1, stats::quantile, c(0.025, 0.5, 0.975)) * 5
         # inc.names = parse.MatrixNames(colnames(ci),numeric=TRUE)
 
         plot(data$time[-1], inc.ci[2, ], type = "n", 
              ylim = range(inc.ci, na.rm = TRUE), ylab = "Ring Increment (mm)")
         PEcAn.visualization::ciEnvelope(data$time[-1], inc.ci[1, ], inc.ci[3, ], col = "lightBlue")
-        points(data$time, data$y[i, ] * 5, pch = "+", cex = 1.5, type = "b", lty = 2)
+        graphics::points(data$time, data$y[i, ] * 5, pch = "+", cex = 1.5, type = "b", lty = 2)
       }
     }
   }
 
   if (FALSE) {
     ## check a DBH
     plot(out[, which(colnames(out) == "x[3,31]")])
-    abline(h = z[3, 31], col = 2, lwd = 2)
-    hist(out[, which(colnames(out) == "x[3,31]")])
-    abline(v = z[3, 31], col = 2, lwd = 2)
+    graphics::abline(h = z[3, 31], col = 2, lwd = 2)
+    graphics::hist(out[, which(colnames(out) == "x[3,31]")])
+    graphics::abline(v = z[3, 31], col = 2, lwd = 2)
   }
 
   ## process model
@@ -59,105 +59,105 @@ InventoryGrowthFusionDiagnostics <- function(jags.out, combined=NULL) {
                            grep("alpha",colnames(out)),
                            grep("deviance",colnames(out)))]
 
-  par(mfrow = c(1, 1))
+  graphics::par(mfrow = c(1, 1))
   for (i in vars) {
-    hist(out[, i], main = colnames(out)[i])
-    abline(v=0,lwd=3)
+    graphics::hist(out[, i], main = colnames(out)[i])
+    graphics::abline(v=0,lwd=3)
   }
   if (length(vars) > 1 && length(vars) < 10) {
-    pairs(out[, vars])
+    graphics::pairs(out[, vars])
   }
 
   if("deviance" %in% colnames(out)){
-    hist(out[,"deviance"])
+    graphics::hist(out[,"deviance"])
     vars <- c(vars,which(colnames(out)=="deviance"))
   }
 
 
   ## rebuild coda for just vars
-  var.out <- as.mcmc.list(lapply(jags.out,function(x){ x[,vars]}))
+  var.out <- coda::as.mcmc.list(lapply(jags.out,function(x){ x[,vars]}))
 
   ## convergence
-  gelman.diag(var.out)
+  coda::gelman.diag(var.out)
 
   #### Diagnostic plots
   plot(var.out)
 
   if("deviance" %in% colnames(out)){
-    hist(out[,"deviance"])
+    graphics::hist(out[,"deviance"])
     vars <- c(vars,which(colnames(out)=="deviance"))
   }
 
   ## rebuild coda for just vars
-  var.out <- as.mcmc.list(lapply(jags.out,function(x){ x[,vars]}))
+  var.out <- coda::as.mcmc.list(lapply(jags.out,function(x){ x[,vars]}))
 
   ## convergence
-  gelman.diag(var.out)
+  coda::gelman.diag(var.out)
 
   #### Diagnostic plots
   plot(var.out)
 
 
   ## Standard Deviations layout(matrix(c(1,2,3,3),2,2,byrow=TRUE))
-  par(mfrow = c(2, 3))
+  graphics::par(mfrow = c(2, 3))
   prec <- out[, grep("tau", colnames(out))]
   for (i in seq_along(colnames(prec))) {
-    hist(1 / sqrt(prec[, i]), main = colnames(prec)[i])
+    graphics::hist(1 / sqrt(prec[, i]), main = colnames(prec)[i])
   }
-  cor(prec)
+  stats::cor(prec)
   # pairs(prec)
 
   ### alpha
-  par(mfrow = c(1, 1))
+  graphics::par(mfrow = c(1, 1))
   alpha.cols <- grep("alpha", colnames(out))
   if (length(alpha.cols) > 0) {
     alpha.ord <- 1:length(alpha.cols)
-    ci.alpha <- apply(out[, alpha.cols], 2, quantile, c(0.025, 0.5, 0.975))
+    ci.alpha <- apply(out[, alpha.cols], 2, stats::quantile, c(0.025, 0.5, 0.975))
     plot(alpha.ord, ci.alpha[2, ], type = "n", 
          ylim = range(ci.alpha, na.rm = TRUE), ylab = "Random Effects")
     PEcAn.visualization::ciEnvelope(alpha.ord, ci.alpha[1, ], ci.alpha[3, ], col = "lightBlue")
-    lines(alpha.ord, ci.alpha[2, ], lty = 1, lwd = 2)
-    abline(h = 0, lty = 2)
+    graphics::lines(alpha.ord, ci.alpha[2, ], lty = 1, lwd = 2)
+    graphics::abline(h = 0, lty = 2)
   }
 
-  par(mfrow = c(1, 1))
+  graphics::par(mfrow = c(1, 1))
   ### alpha
   alpha.cols <- grep("alpha", colnames(out))
   if (length(alpha.cols) > 0) {
     alpha.ord <- 1:length(alpha.cols)
-    ci.alpha <- apply(out[, alpha.cols], 2, quantile, c(0.025, 0.5, 0.975))
+    ci.alpha <- apply(out[, alpha.cols], 2, stats::quantile, c(0.025, 0.5, 0.975))
     plot(alpha.ord, ci.alpha[2, ], type = "n", 
          ylim = range(ci.alpha, na.rm = TRUE), ylab = "Random Effects")
     PEcAn.visualization::ciEnvelope(alpha.ord, ci.alpha[1, ], ci.alpha[3, ], col = "lightBlue")
-    lines(alpha.ord, ci.alpha[2, ], lty = 1, lwd = 2)
-    abline(h = 0, lty = 2)
+    graphics::lines(alpha.ord, ci.alpha[2, ], lty = 1, lwd = 2)
+    graphics::abline(h = 0, lty = 2)
   }
 
   ### YEAR
   year.cols <- grep("year", colnames(out))
   if (length(year.cols > 0)) {
-    ci.yr <- apply(out[, year.cols], 2, quantile, c(0.025, 0.5, 0.975))
+    ci.yr <- apply(out[, year.cols], 2, stats::quantile, c(0.025, 0.5, 0.975))
     plot(data$time, ci.yr[2, ], type = "n", 
          ylim = range(ci.yr, na.rm = TRUE), ylab = "Year Effect")
     PEcAn.visualization::ciEnvelope(data$time, ci.yr[1, ], ci.yr[3, ], col = "lightBlue")
-    lines(data$time, ci.yr[2, ], lty = 1, lwd = 2)
-    abline(h = 0, lty = 2)
+    graphics::lines(data$time, ci.yr[2, ], lty = 1, lwd = 2)
+    graphics::abline(h = 0, lty = 2)
   }
 
   ### INDIV
   ind.cols <- which(substr(colnames(out), 1, 3) == "ind")
   if (length(ind.cols) > 0 & !is.null(combined)) {
-    boxplot(out[, ind.cols], horizontal = TRUE, outline = FALSE, col = as.factor(combined$PLOT))
-    abline(v = 0, lty = 2)
+    graphics::boxplot(out[, ind.cols], horizontal = TRUE, outline = FALSE, col = as.factor(combined$PLOT))
+    graphics::abline(v = 0, lty = 2)
     tapply(apply(out[, ind.cols], 2, mean), combined$PLOT, mean)
     table(combined$PLOT)
 
     spp <- combined$SPP
     # boxplot(out[order(spp),ind.cols],horizontal=TRUE,outline=FALSE,col=spp[order(spp)])
-    boxplot(out[, ind.cols], horizontal = TRUE, outline = FALSE, col = spp)
-    abline(v = 0, lty = 2)
+    graphics::boxplot(out[, ind.cols], horizontal = TRUE, outline = FALSE, col = spp)
+    graphics::abline(v = 0, lty = 2)
     spp.code <- levels(spp)[table(spp) > 0]
-    legend("bottomright", legend = rev(spp.code), col = rev(which(table(spp) > 0)), lwd = 4)
+    graphics::legend("bottomright", legend = rev(spp.code), col = rev(which(table(spp) > 0)), lwd = 4)
     tapply(apply(out[, ind.cols], 2, mean), combined$SPP, mean)
   }
 } # InventoryGrowthFusionDiagnostics

modules/data.land/R/extract_FIA.R

@@ -20,7 +20,7 @@ extract_FIA <- function(lon, lat, start_date, end_date, gridres = 0.075, dbparms
   veg_info <- list()
 
   fia.con <- PEcAn.DB::db.open(dbparms$fia)
-  on.exit(db.close(fia.con), add = TRUE)
+  on.exit(PEcAn.DB::db.close(fia.con), add = TRUE)
 
   lonmin   <- lon - gridres
   lonmax   <- lon + gridres
@@ -49,7 +49,7 @@ extract_FIA <- function(lon, lat, start_date, end_date, gridres = 0.075, dbparms
                    " AND p.lat <= ", latmax, " AND p.measyear >= ", min.year, 
                    " AND p.measyear <= ", max.year, " GROUP BY p.cn")
 
-  plot.info <- db.query(query, con = fia.con)
+  plot.info <- PEcAn.DB::db.query(query, con = fia.con)
   if (nrow(plot.info) == 0) {
     PEcAn.logger::logger.severe("No plot data found on FIA.")
   }
@@ -97,7 +97,7 @@ extract_FIA <- function(lon, lat, start_date, end_date, gridres = 0.075, dbparms
                     " and p.lon < ", lonmax, 
                     " and p.lat >= ", latmin,
                     " and p.lat < ", latmax)
-  tree.info <- db.query(query, con = fia.con)
+  tree.info <- PEcAn.DB::db.query(query, con = fia.con)
   names(tree.info) <- tolower(names(tree.info))
 
   if (nrow(tree.info) == 0) {