When gdm.varImp() drops a predictor variable due to zero sum of I-spline coefficients, it throws a dimname error when setting the colnames of varImpTable

Error:

GPC.varImp <- gdm::gdm.varImp(GPC_dat, geo = T, nPerm = 1000, parallel = FALSE)
Fitting initial model with all 4 predictors...
Sum of I-spline coefficients for predictors(s) X = 0
Removing predictor(s) X and proceeding with permutation testing...
Error in dimnames(x) <- dn : 
  length of 'dimnames' [2] not equal to array extent

Suggested fix is to change lines 411-416 of gdm.variable.importance.R from:

  varImpTable <- matrix(NA, nVars, nVars*5-1)
  rownames(varImpTable) <- varNames
  colnames(varImpTable) <- c("All predictors",
                             "1-removed",
                             paste(seq(2,nVars*5-2), "-removed", sep=""))

to something like:

  varImpTable <- matrix(NA, nVars, nVars - 1)
  rownames(varImpTable) <- varNames
  varImpTable_colnames <- c("All predictors", "1-removed")
  if(nVars - 2 > 1) {
    varImpTable_colnames <- c(varImpTable_colnames, paste(seq(2, nVars - 2), "-removed", sep = ""))
  }
  colnames(varImpTable) <- varImpTable_colnames

sessioninfo::session_info()
#> ─ Session info ───────────────────────────────────────────────────────────────
#>  setting  value
#>  version  R version 4.2.0 (2022-04-22)
#>  os       macOS Monterey 12.6
#>  system   aarch64, darwin20
#>  ui       X11
#>  language (EN)
#>  collate  en_US.UTF-8
#>  ctype    en_US.UTF-8
#>  tz       America/Los_Angeles
#>  date     2022-09-27
#>  pandoc   2.18 @ /Applications/RStudio.app/Contents/MacOS/quarto/bin/tools/ (via rmarkdown)
#> 
#> ─ Packages ───────────────────────────────────────────────────────────────────
#>  package     * version date (UTC) lib source
#>  cli           3.3.0   2022-04-25 [1] CRAN (R 4.2.0)
#>  clipr         0.8.0   2022-02-22 [1] CRAN (R 4.2.0)
#>  digest        0.6.29  2021-12-01 [1] CRAN (R 4.2.0)
#>  evaluate      0.16    2022-08-09 [1] CRAN (R 4.2.0)
#>  fastmap       1.1.0   2021-01-25 [1] CRAN (R 4.2.0)
#>  fs            1.5.2   2021-12-08 [1] CRAN (R 4.2.0)
#>  glue          1.6.2   2022-02-24 [1] CRAN (R 4.2.0)
#>  highr         0.9     2021-04-16 [1] CRAN (R 4.2.0)
#>  htmltools     0.5.3   2022-07-18 [1] CRAN (R 4.2.0)
#>  knitr         1.40    2022-08-24 [1] CRAN (R 4.2.0)
#>  lifecycle     1.0.1   2021-09-24 [1] CRAN (R 4.2.0)
#>  magrittr      2.0.3   2022-03-30 [1] CRAN (R 4.2.0)
#>  reprex        2.0.2   2022-08-17 [1] CRAN (R 4.2.0)
#>  rlang         1.0.5   2022-08-31 [1] CRAN (R 4.2.0)
#>  rmarkdown     2.16    2022-08-24 [1] CRAN (R 4.2.0)
#>  rstudioapi    0.14    2022-08-22 [1] CRAN (R 4.2.0)
#>  sessioninfo   1.2.2   2021-12-06 [1] CRAN (R 4.2.0)
#>  stringi       1.7.8   2022-07-11 [1] CRAN (R 4.2.0)
#>  stringr       1.4.1   2022-08-20 [1] CRAN (R 4.2.0)
#>  withr         2.5.0   2022-03-03 [1] CRAN (R 4.2.0)
#>  xfun          0.32    2022-08-10 [1] CRAN (R 4.2.0)
#>  yaml          2.3.5   2022-02-21 [1] CRAN (R 4.2.0)
#> 
#>  [1] /Library/Frameworks/R.framework/Versions/4.2-arm64/Resources/library
#> 
#> ──────────────────────────────────────────────────────────────────────────────

I run gdm.varImp() with only one covariate (I tried both "awcA" and "phTotal") and it returns "NA" for "model deviance", "p-value" etc. This does not happen if I run it with both "awcA" and "phTotal".

library(gdm)

sppData <- southwest[c(1,2,13,14)]
envTab <- southwest[c(2:ncol(southwest))]
sitePairTab <- formatsitepair(sppData, 2, XColumn="Long", YColumn="Lat",
                              sppColumn="species", siteColumn="site", predData=envTab)

library(tidyverse)
test<- sitePairTab %>% select(distance,weights,s1.xCoord,s1.yCoord,s2.xCoord,
                              s2.yCoord,s1.awcA,s2.awcA)


modTest <- gdm.varImp(test, geo=F, nPerm=50, parallel=T, cores=10)

> modTest$`Predictor p-values`
     All predictors
awcA             NA
> modTest$`Model assessment`
                           All predictors
Model deviance                         NA
Percent deviance explained             NA
Model p-value                          NA
Fitted permutations                    NA
> modTest$`Predictor Importance`
     All predictors
awcA             NA
> modTest$`Model Convergence`
     All predictors
awcA             NA
> modTest$`Predictor p-values`
     All predictors
awcA             NA
> modTest$`Model assessment`
                           All predictors
Model deviance                         NA
Percent deviance explained             NA
Model p-value                          NA
Fitted permutations                    NA
> modTest$`Predictor Importance`
     All predictors
awcA             NA
> modTest$`Model Convergence`
     All predictors
awcA             NA

The original isplineExtract() function is hardcoded for the default of 3 splines per variable. The problem is in this line of the function definition:

xDat[,i] <-  seq(from=model$knots[[(i*3)-2]],to=model$knots[[(i*3)]], length=PSAMPLE)

And it can be fixed with:

xDat[,i] <-  seq(from=model$knots[[(i*numsplines)-(numsplines-1)]],to=model$knots[[(i*numsplines)]], length=PSAMPLE)

I followed the same code used in vignette (https://github.com/fitzLab-AL/GDM/blob/master/vignettes/gdmVignette.pdf) for formatsiterepair(). Although my data structure was exactly similar to that used in the vignette. My dataset has about 6000 rows. Using 1000 rows don't give any error. There is not NA values in my dataset. Any solution please?

> gdmTab <- formatsitepair(sppTab, bioFormat = 2, siteColumn = "site", XColumn = "Long", YColumn = "Lat", sppColumn = "grasshopper", predData = envTab) Aggregation function missing: defaulting to length Error in createsitepair(dist = distData, spdata = bioData, envInfo = predData, : The length of distance values are not the same as the expected number of rows of the site-pair table, unable to proceed. In addition: Warning message: In formatsitepair(sppTab, bioFormat = 2, siteColumn = "site", XColumn = "Long", : No abundance column specified, assuming species data are presence

It is not a issue but I would like to open a discussion on one assumption I found in the Ferrier et al (2007), they say:

In GDM we make the reasonable assumption that compositional dissimilarity can only increase, not decrease, with increasing separation of sites along an environmental gradient.

Despite I think it is a reasonable assumption I was wondering that the opposite may be easy to test if we scale the data to its inverse. For instance, first I can check if plant dissimilarity increases in a increasing P-content in soil and then I would look if it increase in an scale that reverses the pH values. Does it make sense?

Another approach would be to analyze the data in a similarity distance instead of a dissimilarity one. What would be the best approach?

From line 286 in gdm.variable.importance.R predictors got removed where the I-Spline coefficients == 0. If this is the "Geographic" predictor which is not a column name in the spTable, the whole table gets empty. I don't have the full overview for this rather complex function, but I think it would be enough to check for this case and set geo <- F:

  zeroSum <- fullGDM$predictors[which(sumCoeff==0)]
  if(length(zeroSum)>0){
    message(paste0("Sum of I-spline coefficients for predictors(s) ", zeroSum," = 0"))
    Sys.sleep(1)
    message(paste0("Removing predictor(s) ", zeroSum, " and proceeding with permutation testing..."))
    Sys.sleep(1)
    
    for(z in zeroSum){
      ##select variable columns to be removed from original site-pair table
      if(z == "Geographic") { geo <- F} else {
        testVarCols1 <- grep(paste("^s1.", z, "$", sep=""), colnames(spTable))
        testVarCols2 <- grep(paste("^s2.", z, "$", sep=""), colnames(spTable))
        spTable <- spTable[,-c(testVarCols1, testVarCols2)]
      }
    }
  }

Hi @fitzLab-AL @mliskAL

I would like to make some predictions but got an error.
I used two sets of factors in my gdm model,1) 19 bioclimate factors(original value) and 2) land cover data(range from 0~1). When I tried to make a prediction through time the result only had zero value, I checked whether this was brought out by having the same raster layer(i.e. climate or land cover not change) for the two time periods stack but the answer was not.
I've had no problems with the gdm model using 19 bioclimate factors in the past. So I run two models with the same response data but using only climate or land cover data, this time the predict function worked well.
So I cannot figure out why the function produces an error result when the two sets of factors are combined, they have the same resolution and extent, and different period stacks have the same order for factors.

Looking forward to your reply and thanks in advance.

Hello @fitzLab-AL and gdm users (also @rvalavi if you're inclined),

Is the code for reproducing plot 4a available? Section 1.7 of the supplementary materials (and online material - link below) includes code for reproducing plot 4b, but not 4a. A similar plot (to 4a) is one of many produced via plot(gdm). I can reproduce that individual plot from plot(gdm) (in part) using plot(gdm$ecological, gdm$observed). But I lack the line showing GDM-predicted dissimilarity.

Having this code will also help clarify whether a data object from predict(gdm, gdmTab) is required to produce plot 4a; and thus further reveal the relationship between plots produced via plot(gdm) and plot(x, gdm.pred).

My expectation is that other users have had this question. Thanks very much

[https://github.com/fitzLab-AL/gdm?tab=readme-ov-file#visualizing-multi-dimensional-biological-patterns]

Due to line 569 in the file gdm.variable.importance.R:
grepper <- grep(var, names(permVarDev))
the function gdm.varImp failes if the varible names are character subssets of each other. For example for varibles "var1" and "var10" grep returns both while searching for the first one. The function failes afterwards while trying:
varDevTab <- permVarDev[[grepper]]
with the error message: "Error in permVarDev[[grepper]] : recursive indexing failed at level 2"
I would suggest using:
grepper <- which(var == names(permVarDev))
at line 569 to search for an exact match.

I'm trying to use formatsitepair to analyse the effect of two sets of distance matrices on a third one, but the formatsitepair function gives the following error:
Error in formatsitepair(bioData = In.fmt, bioFormat = 3, abundance = TRUE, :
argument "predData" is missing, with no default

However, in the help of this function, distPreds is stated to be:
An optional list of distance matrices to be used as predictors, either in combination with, or as a substitute for, predData.

How am I supposed to use only distance matrices and no predData?

I have found that the results of a GDM can be very sensitive to spline number (supplied in the splines argument), and I would suggest that future versions of gdm() require users to supply the splines argument to avoid false confidence in the default settings.

Along the same lines, it would be helpful to include some means of comparing models of different complexity to decide what the spline number should be. Ferrier et al. (2007) suggest such procedures, albeit vaguely:

"The predictors included in a model, and the number of I-spline basis functions employed for each predictor (controlling the allowable complexity of the nonlinear transformation of that predictor), can be determined using various automated selection strategies, including forward-selection, backward-elimination, and stepwise procedures."

Hello, I would like to perform GDM, but my data comes from a designed experiment so the geographical distance do not provide much information. Therefore, I would like to know if it is possible to apply GDM without adding this information.

Thanks in advance

It took me a while to figure out that the plots of predictor variable splines with plot.gdm() were done sequentially, and that I had to press BACK in RStudio to find them. A more intuitive way would be if plot.gdm() set par(ask=TRUE) before plotting to allow users to click through the plots, then par(ask=FALSE) to turn it off at the end.

Dear Fitzpatrick Lab,

I noticed that the geographic distance between two lat-lon points is calculated using the Pythagorean theorem (e.g. v <- sqrt((data[,3]-data[,5])^2 + (data[,4]-data[,6])^2)) rather than a formula that takes into account the spherical(-ish) shape of the Earth. This causes the calculated distances to become more distorted the farther samples are from the equator. May I recommend using something like the haversine formula to calculate geographic distance as an alternative?

Cheers,
Clara

regardless of what value is provided for the subsample argument, function always returns 100% of site pairs.

Hi! I'm trying to produce a gdm with presence/absence data that I have and this is the issue that I am having.

This is the command that I am running, and the error that I get:
gdmtab <- formatsitepair(bioData = spp, bioFormat = 2, XColumn="Longitude", YColumn="Latitude", sppColumn = "Trait", site= "zone", predData = env) No abundance column was specified, so the species data are assumed to be presences.Aggregation function missing: defaulting to length cannot xtfrm data framesError in[.data.frame(envInfo, , "siteUltimateCoolness") : undefined columns selected
Here is the head of my biodata:
head(spp)
Longitude Latitude Trait zone
1 -110.3976 31.46074 n A
2 -104.0201 30.65312 n W
3 -109.2110 31.89133 y B
4 -109.2393 31.87379 y B
5 -107.4298 28.77738 y BB
6 -105.7348 28.92134 n CC

Trait is the y/n of having the trait, and the zone is the locality found. My environment data is worldclim, elevation, Latitude, Longitude. Do you have any suggestions for how to address this issue? Thanks!

Running gdm.varImp() with the example works as expected

library(gdm)

sppData <- southwest[c(1,2,13,14)]
envTab <- southwest[c(2:ncol(southwest))]
sitePairTab <- formatsitepair(sppData, 2, XColumn="Long", YColumn="Lat",
                              sppColumn="species", siteColumn="site", predData=envTab)

modTest <- gdm.varImp(sitePairTab, geo=T, nPerm=50, parallel=T, cores=10)

However, when I try to run it with only a few covariates, it returns an error:

library(tidyverse)
test<- sitePairTab %>% select(distance,weights,s1.xCoord,s1.yCoord,s2.xCoord,
                              s2.yCoord,s1.awcA,s1.phTotal,s2.awcA,s2.phTotal)


modTest <- gdm.varImp(test, geo=F, nPerm=50, parallel=T, cores=10)



Fitting initial model with all 2 predictors...
Error in matrix(NA, 4, nVars, dimnames = list(c("Model deviance", "Percent deviance explained",  : 
  length of 'dimnames' [2] not equal to array extent

Am I doing something wrong?

Hi all,

I would like to assess variance importance using gdm.varImp. In the manual, it says that if user supplied knots are provided, the knots argument should be the same as used in model fitting. However, this results in an error:

EXPLO_gdmRes <- gdm(data = EXPLO_spt,
                          geo = FALSE,
                          splines = EXPLO_splines_vector,
                          knots = EXPLO_knots_values)
      
      # OPERATIONS TO PERFORM IF A MODEL WAS SUCCESSFULLY FIT
      if (!is.null(EXPLO_gdmRes)) {
        # 1. PERFORM PREDICTOR AND MODEL SIGNIFICANCE TESTING
        EXPLO_varImp <- gdm.varImp(
          EXPLO_spt,
          geo = FALSE,
          splines = EXPLO_splines_vector,
          knots = EXPLO_knots_values,
          predSelect = TRUE,
          nPerm = 100,
          pValue = 0.05,
          parallel = TRUE,
          cores = 3,
          sampleSites = 1, # the % of sites to use for purmutation routines; 1 = 100 %
          sampleSitePairs = 1, # the % of site pairs to use for purmutation routines; 1 = 100 %
          outFile = paste0(i, "_", y, "_", s, "_varImp_log.RData")
          )

Creating 100 permuted site-pair tables...
Starting model assessment...
Error in gdm(currSitePair, geo = geo, splines = splines, knots = knots) : 
  When knots are supplied by the user, there should be 6 items in the knots argument, not 141 items.

Obviously, following the construction of permuted site-pair tables, the knots argument is not longer the one used for model fitting.
How can I solve this issue?

Thank your for your assistance.

I would like to use the custom weights option of the formatsitepair function. However the code produces the error: Error in custWeights[, "gettingCoolSiteColumn"][, 1] : incorrect number of dimensions. I think the offending line is hwap <- custWeights[, "gettingCoolSiteColumn"][, 1], which seems to incorrectly subset the custom weights data.frame (unless I've incorrectly set up the custom weights data.frame). Here's a reproducible example:

# set up data
sppData <- southwest[c(1,2,13,14)]
envTab <- southwest[c(2:ncol(southwest))]

# create custom weights data.frame
custw <- data.frame(site = unique(sppData$site), weights = runif(length(unique(sppData$site))))

# run format site pair
sitePairTab <- formatsitepair(sppData, 2, XColumn = "Long", YColumn = "Lat", sppColumn = "species", 
siteColumn = "site", predData = envTab, weightType = "custom", custWeights = custw, verbose = TRUE)

Best,
Rob

I'm sure you're aware of this, but CRAN doesn't allow installing the GDM package anymore. It'd be nice if you updated the package so that it continues to be supported.

Hello,

I have a large gridded (raster) geographic area (>17000 sites) for which I have a phylogenetic distance matrix (phyloSimpson). When I create the site x pair matrix with the distance matrix, the lat/lon info, and my predictor variables, the total vector length of the matrix is over 2^31 elements, and thus is long vector. When I run gdm on this site x pair matrix, I get the following error: "long vectors (argument 2) are not supported in .C" (argument 2 would be the matrix input for GDM_FitFromTable).

I notice the R code for the various functions in this package use .C to call the C++ code. From searching around, it looks like .C and .Fortran cannot handle long vectors but .Call can. When I make a duplicate function locally of gdm using .Call instead of .C (for GDM_FitFromTable), I am able to run the code on my long vector, however, I run into a segmentation fault where memory does not map. I'm guessing this is related to the mismatch in my R environment and the original .cpp from the package, though I am not sure.

Would it be possible to make an adjustment for long vectors in the gdm package? I could also subset the matrix and use predict.gdm for the other cells if needed, though I would like to try gdm for the entire matrix, and long vectors may be useful for other large raster-based projects.

Please let me know if you would like any more information or if I am mistaken about the source of the issue. Thank you.

Hello,
gdm.varImp() is returning the following error message:

"Fitting initial model with all 18 predictors...
Sum of I-spline coefficients for predictors(s) Geographic = 0Sum of I-spline coefficients for predictors(s) t.mean = 0Sum of I-spline coefficients for predictors(s) t.max = 0
Removing predictor(s) Geographic and proceeding with permutation testing...Removing predictor(s) t.mean and proceeding with permutation testing...Removing predictor(s) t.max and proceeding with permutation testing...
Error in [.data.frame(spTable, c(7:(6 + nVars))) :
undefined columns selected"

Used this function before with no issues. Any hints? Thanks a lot!

Hello,
I've encountered two issues with the gdm.varImp() function. The first is easy to solve:
For GDMs fit to just one environmental variable, with geo=T, gdm.varImp() fails to run, due to the line

varNames <- colnames(spTable[,c(7:6+nVars))])

Presumably this should be replaced with

varNames <- colnames(spTable)[c(7:6+nVars))]

Second:
gdm.varImp() decides how many sites to expect based on the unique number of site-coordinates. This seems reasonable enough, except that some users (like me) will input two different site names with different environmental conditions at a single set of coordinates (due to a matched-pair design). This is easy enough to solve by adding a tiny offset to the coordinates of the sites in condition 2, but it would probably stump anybody not comfortable looking at the source code. It would be nice if gdm.varImp figured out how many sites to expect based on the site names, instead.

Cheers
Jacob

A couple notes on plotUncertainty():

1. After reading the documentation, I'm still a bit unclear on exactly how exactly it works and how I would explain it in a paper. Could the documentation on this command be elaborated a bit?
2. I use isplineExtract() to tabulate splines for custom plotting, but there does not seem to be any corresponding method for extract the results of plotUncertainty() for custom plotting. This would be nice.

This works fine:
gdm.varImp_WN <- gdm.varImp(sitepair_WN, geo = FALSE)

This breaks in a puzzling way.
gdm.varImp_WN <- gdm.varImp(sitepair_WN, geo = FALSE, splines = c(6, 6, 6, 6, 6))

Error in gdm(currSitePair, geo = geo, splines = splines, knots = knots) : 
  Number of splines does not equal the number of predictors. 
              Splines argument has 6 items but needs 5 items.

One might think that shortening the splines argument by one item would satisfy it:
gdm.varImp_WN <- gdm.varImp(sitepair_WN, geo = FALSE, splines = c(6, 6, 6, 6))

. . . but it is fickle. Now it wants the 5th item back.

Error in gdm(currSitePair, geo = geo, splines = splines, knots = knots) : 
  Number of splines does not equal the number of predictors. 
              Splines argument has 4 items but needs 5 items.

R-Studio crashes when trying to run the gdm.transform function for a rasterStack object, giving the following error:

Error in transformed[, i] <- tmp :
'pairlist' object cannot be coerced to type 'double'

Hello,

I am using the latest CRAN version of GDM (v. 1.4.2.2) and I am able to run the gdm function on my data and the results look reasonable.

However, I keep getting this error when running gdm.varImp:

modTest_b <- gdm.varImp(gdm_sitepair, fullModelOnly = F, geo=T, nPerm=100, parallel=F,splines=rep(3,13))
Error in data[, 3] - data[, 5] : non-numeric argument to binary operator

I realize that asking for a solution without details of the input data is a bit nebulous.

Therefore, I have attached the data that can be read in with dget(gdm_sitepair.txt)->gdm.varImp.txt

Advice would be really appreciated!

Thanks so much for any help!

Frank

gdm_sitepair.txt

Hi,
I'm running Ubuntu 20.04 and R 4.1.3 and I'm getting this installation error. I also tried devtools::install_github("fitzLab-AL/gdm") with the same result.

install.packages("gdm")
Installing package into ‘/home/fd/R/x86_64-pc-linux-gnu-library/4.1’
(as ‘lib’ is unspecified)
trying URL 'https://cloud.r-project.org/src/contrib/gdm_1.5.0-1.tar.gz'
Content type 'application/x-gzip' length 688341 bytes (672 KB)
==================================================
downloaded 672 KB

* installing *source* package ‘gdm’ ...
** package ‘gdm’ successfully unpacked and MD5 sums checked
** using staged installation
** libs
g++ -std=gnu++11 -I"/usr/share/R/include" -DNDEBUG  -I'/usr/lib/R/site-library/Rcpp/include'    -fpic  -g -O2 -fdebug-prefix-map=/build/r-base-lENDSu/r-base-4.1.3=. -fstack-protector-strong -Wformat -Werror=format-security -Wdate-time -D_FORTIFY_SOURCE=2 -g  -c Gdmlib.cpp -o Gdmlib.o
g++ -std=gnu++11 -I"/usr/share/R/include" -DNDEBUG  -I'/usr/lib/R/site-library/Rcpp/include'    -fpic  -g -O2 -fdebug-prefix-map=/build/r-base-lENDSu/r-base-4.1.3=. -fstack-protector-strong -Wformat -Werror=format-security -Wdate-time -D_FORTIFY_SOURCE=2 -g  -c Message.cpp -o Message.o
g++ -std=gnu++11 -I"/usr/share/R/include" -DNDEBUG  -I'/usr/lib/R/site-library/Rcpp/include'    -fpic  -g -O2 -fdebug-prefix-map=/build/r-base-lENDSu/r-base-4.1.3=. -fstack-protector-strong -Wformat -Werror=format-security -Wdate-time -D_FORTIFY_SOURCE=2 -g  -c NNLSDoubleCore.cpp -o NNLSDoubleCore.o
g++ -std=gnu++11 -I"/usr/share/R/include" -DNDEBUG  -I'/usr/lib/R/site-library/Rcpp/include'    -fpic  -g -O2 -fdebug-prefix-map=/build/r-base-lENDSu/r-base-4.1.3=. -fstack-protector-strong -Wformat -Werror=format-security -Wdate-time -D_FORTIFY_SOURCE=2 -g  -c NNLS_Double.cpp -o NNLS_Double.o
g++ -std=gnu++11 -I"/usr/share/R/include" -DNDEBUG  -I'/usr/lib/R/site-library/Rcpp/include'    -fpic  -g -O2 -fdebug-prefix-map=/build/r-base-lENDSu/r-base-4.1.3=. -fstack-protector-strong -Wformat -Werror=format-security -Wdate-time -D_FORTIFY_SOURCE=2 -g  -c dllmain.cpp -o dllmain.o
gcc -I"/usr/share/R/include" -DNDEBUG  -I'/usr/lib/R/site-library/Rcpp/include'    -fpic  -g -O2 -fdebug-prefix-map=/build/r-base-lENDSu/r-base-4.1.3=. -fstack-protector-strong -Wformat -Werror=format-security -Wdate-time -D_FORTIFY_SOURCE=2 -g  -c registerDynamicSymbol.c -o registerDynamicSymbol.o
g++ -std=gnu++11 -I"/usr/share/R/include" -DNDEBUG  -I'/usr/lib/R/site-library/Rcpp/include'    -fpic  -g -O2 -fdebug-prefix-map=/build/r-base-lENDSu/r-base-4.1.3=. -fstack-protector-strong -Wformat -Werror=format-security -Wdate-time -D_FORTIFY_SOURCE=2 -g  -c stdafx.cpp -o stdafx.o
g++ -std=gnu++11 -shared -L/usr/lib/R/lib -Wl,-Bsymbolic-functions -Wl,-z,relro -o gdm.so Gdmlib.o Message.o NNLSDoubleCore.o NNLS_Double.o dllmain.o registerDynamicSymbol.o stdafx.o -L/usr/lib/R/lib -lR
installing to /home/fd/R/x86_64-pc-linux-gnu-library/4.1/00LOCK-gdm/00new/gdm/libs
** R
** data
*** moving datasets to lazyload DB
** inst
** byte-compile and prepare package for lazy loading
double free or corruption (out)
Aborted (core dumped)
ERROR: lazy loading failed for package ‘gdm’
* removing ‘/home/fd/R/x86_64-pc-linux-gnu-library/4.1/gdm’
Warning in install.packages :
  installation of package ‘gdm’ had non-zero exit status

The downloaded source packages are in
	‘/tmp/RtmpfbOcEj/downloaded_packages’

I would like two add two columns to an object created by formatsitepair() containing unique identifiers for the sites that are being compared by the Bray-Curtis index. Something like this:

site_id1, site_id2, distance, weights,
1,2,0.5,1
1,3,0.8,1

I run the example code

head(southwest)
sppData <- southwest[c(1,2,13,14)]
envTab <- southwest[c(2:ncol(southwest))]

sitePairTab <- formatsitepair(sppData, 2, XColumn="Long", YColumn="Lat", sppColumn="species",
                              siteColumn="site", predData=envTab)

head(sitePairTab )

 distance weights s1.xCoord s1.yCoord s2.xCoord s2.yCoord s1.awcA s1.phTotal s1.sandA  s1.shcA s1.solumDepth
132   0.4485981       1   115.057 -29.40472  115.5677 -29.46599 23.0101   480.3266 83.99326 477.5656      1129.933
132.1 0.7575758       1   115.057 -29.40472  116.0789 -29.52556 23.0101   480.3266 83.99326 477.5656      1129.933
132.2 0.8939394       1   115.057 -29.40472  116.5907 -29.58342 23.0101   480.3266 83.99326 477.5656      1129.933
132.3 0.9178082       1   115.057 -29.40472  117.1029 -29.63957 23.0101   480.3266 83.99326 477.5656      1129.933

and noticed that some row.names seem to be a combination of two numbers (e.g. 131.1) but others do not (e.g. 132). I also noticed these numbers do not to match the column "site" from "envTab".

I think I need help from someone who knows the inner workings of the function.

Thank you.

When I attempt to run gdm.varimp() I am returned with an error "Error in [.data.frame(x, , -c(remVarCols1, remVarCols2)) :
object 'remVarCols1' not found"

I have a number of datasets and for some of them the code works and for others it doesn't. The only difference is the species included in the analysis, the landscape variables are all the same for each site. I have attached an example to demonstrate the error, I hope it works! (and by that I mean I hope it gives you an error)

Thank you

Hello @rvalavi and @fitzLab-AL,

I'm testing gdm ver 1.6.0 using terra data objects. To lower computational time for this test, I limited my model to 3 predictors (instead of the 17 I'm normally running) with geo=TRUE. Tests were run on a desktop PC.

The code in Section 2.1 from Appendix S1 of the user guide (Mokany et al 2022) ran okay with a couple minor code amendments. Next, I want to create maps so I can visualize compositional variation across geographic space. The code in section 2.2 of Appendix S1 worked fine until I got to "scaling the PCA rasters to make full use of the colour spectrum"

I get an error when running any one of the following code chunks. Note, I added a space after the "at" symbol in these code chunks because a github user has that as a handle i.e., @ data (no space)

pcaRast[[1]] <- (pcaRast[[1]]-pcaRast[[1]]@ data@min) / (pcaRast[[1]]@ data@max-pcaRast[[1]]@ data@min)*255
pcaRast[[2]] <- (pcaRast[[2]]-pcaRast[[2]]@ data@min) / (pcaRast[[2]]@ data@max-pcaRast[[2]]@ data@min)*255
pcaRast[[3]] <- (pcaRast[[3]]-pcaRast[[3]]@ data@min) / (pcaRast[[3]]@ data@max-pcaRast[[3]]@ data@min)*255

-- Error: no slot of name "data" for this object of class "SpatRaster"

I'm assuming this is because SpatRaster objects were not part of the original pipeline. I tried changing the references (e.g., [[1]]) in these chunks and it didn't work. Note, the name of each raster in pcaRast are lyr1, lyr2, and lyr3. Any suggestions?

Overall, the addition of terra() functionality to this version of gdm makes things much faster, which is great to see.

Thanks.

Hello all,
first of all, thank you for providing this powerful piece of software to the scientific community.

I have quite complex models to fit for multiple species and sub-setting approaches. If I use the default settings (# splines 3), at least two models can be fit. However, as increasing the number of splines and knots might result in accommodating more complex models, I hoped I could get more models to fit by increasing these parameters. From the manual, it is not completely clear whether the number of knots is automatically adjusted for the number of splines requested. Hence, I decided to calculate the quantiles for the knots myself using the following approach:

 EXPLO_npred_splines <- ((ncol(EXPLO_spt) - 6) / 2)
      EXPLO_number_splines <- GLOBAL_nsplines
      EXPLO_splines_vector <-
        rep(EXPLO_number_splines, EXPLO_npred_splines)
      
      # Adjust knot vector
      EXPLO_knots_values <- c()
      for (o in 2:ncol(EXPLO_env_data)) {
        EXPLO_pred_col <- EXPLO_env_data[o]
        EXPLO_quant_vec <-
          c(0:(EXPLO_number_splines - 1)) / (EXPLO_number_splines - 1)
        EXPLO_pred_vals <-
          quantile(EXPLO_pred_col,
                   EXPLO_quant_vec,
                   na.rm = TRUE,
                   names = FALSE)
        # Append the knots values to vector to hold the values for all predictors
        EXPLO_knots_values <-
          c(EXPLO_knots_values, EXPLO_pred_vals)
      }

However, if I pursue this approach, the two models that were fit with the default values fail. I tried different values for the number of splines (4, 5, 10, 50, 100) but still, no model can be fit. I am aware of the fact that the impact of these parameters is not fully understood yet (https://github.com/fitzLab-AL/gdm#gdm-fitting). However, I am really surprised that even changing the number of splines to 4 results in model fit failures and I still don't understand why.

Thank you for your assistance.

Dear Fitzpatrick Lab,

First, thanks for this great package. I am experiencing an issue when running the function formatsitepair. As biological data, I use a dissimilarity matrix, and as environmental predictor I use a distance matrix with the argument distPreds. As the formatsitepair function does not accept only distances matrices as predictors, I provide one fake predictor for the argument predData what includes coordinate columns and a fake variable. However, it returns the following error:

Error in if ((dim(mat1)[1] != dim(mat)[1]) & (dim(mat1)[2] != dim(mat)[2])) { :
argument is of length zero

This is my code, and I attach the input files:

sim <- read.table("sim.txt", h=T)
euc= read.table("euc.txt", h=T)
fake= read.table("fake.txt", h=T)

simdata = formatsitepair(bioData=sim, bioFormat=3, predData=fake, distPreds=euc, XColumn="x", YColumn="y", siteColumn="site")

I am not able to find what I am doing wrong here.

Many tanks in advance,

Pedro

euc.txt
fake.txt
sim.txt

When I specify non-equal numbers of splines for the variables in my model, the model call itself seems to run correctly, but both plot() and isplineExtract() throw "subscript out of bounds" errors. I would presume this is because these functions are expecting equal dimensions across variables and then trying to access rows that don't exist.

gdm_WN6 <- gdm(sitepair_WN, geo = FALSE, splines = c(6, 4, 5, 6, 6))

plot(gdm_WN6)

Error in h(simpleError(msg, call)) : 
  error in evaluating the argument 'x' in selecting a method for function 'plot': subscript out of bounds

Hi, @fitzLab-AL & @mliskAL

I have trouble running phylo-GDM, the site-by-site biological distance (dissimilarity) is either 0 or 1, how to address it?
Can gdm package implement cluster(Laffan et al., 2010) based on a Sorenson distance in Biodiverse tool?

Data:
Example_tree.nex
Example_site_data.csv

library(gdm)
library(ape)
library(tidytree)
tre <- read.nexus("Example_tree.nex")
tre$tip.label <- gsub("Genus:","",tre$tip.label)
tre_tb <- as_tibble(tre)
dt <- read.csv("Example_site_data.csv")
test <- dt %>% left_join(tre_tb[,3:4], by = c("species" = "label"))

sp <- test[,c("num","species","x", "y")]
phy <- test[,c("num","branch.length")]

exFormat <- formatsitepair(sp, bioFormat=2, siteColumn="num", XColumn="x",sppColumn = "species",
                             YColumn="y", predData=phy)

head(exFormat)

distance weights s1.xCoord s1.yCoord s2.xCoord s2.yCoord s1.branch.length s2.branch.length
1        0       1   3229628   3078198   3216986   2951192       0.07551968       0.16666667
2        0       1   3229628   3078198   3216038   2960749       0.07551968       0.07551968
3        1       1   3229628   3078198   3217265   2960874       0.07551968       0.16666667
4        1       1   3229628   3078198   3205746   2963060       0.07551968       0.07551968
5        1       1   3229628   3078198   3214891   2928147       0.07551968       0.16666667
6        1       1   3229628   3078198   3219317   2971030       0.07551968       0.03299437

reference
Laffan, S. W., Lubarsky, E., & Rosauer, D. F. (2010). Biodiverse, a tool for the spatial analysis of biological and related diversity. Ecography, 33(4), 643-647. 10.1111/j.1600-0587.2010.06237.x

Thank you！

Looking forward to your reply.