CRAN Package Check Results for Maintainer ‘Alexander Ilich <ailich at usf.edu>’

Last updated on 2025-12-03 05:50:38 CET.

Package	ERROR	NOTE	OK
GLCMTextures	1		12
MultiscaleDTM	2	2	9

Package GLCMTextures

Current CRAN status: ERROR: 1, OK: 12

Version: 0.6.3
Check: examples
Result: ERROR Running examples in 'GLCMTextures-Ex.R' failed The error most likely occurred in: > ### Name: glcm_textures > ### Title: Calculates GLCM texture metrics of a Raster Layer > ### Aliases: glcm_textures > > ### ** Examples > > r<- rast(volcano, extent= ext(2667400, 2667400 + ncol(volcano)*10, + 6478700, 6478700 + nrow(volcano)*10), crs = "EPSG:27200") Warning: PROJ: proj_create_from_database: Cannot find proj.db (GDAL error 1) Error: [rast] empty srs Execution halted Flavor: r-oldrel-windows-x86_64

Version: 0.6.3
Check: tests
Result: ERROR Running 'testthat.R' [13s] Running the tests in 'tests/testthat.R' failed. Complete output: > # This file is part of the standard setup for testthat. > # It is recommended that you do not modify it. > # > # Where should you do additional test configuration? > # Learn more about the roles of various files in: > # * https://r-pkgs.org/testing-design.html#sec-tests-files-overview > # * https://testthat.r-lib.org/articles/special-files.html > > library(testthat) > library(GLCMTextures) Loading required package: terra terra 1.8.86 Attaching package: 'terra' The following objects are masked from 'package:testthat': compare, describe > > test_check("GLCMTextures") Saving _problems/test-glcm_textures-45.R Saving _problems/test-glcm_textures-54.R Saving _problems/test-glcm_textures-63.R Saving _problems/test-glcm_textures-75.R [ FAIL 4 | WARN 4 | SKIP 0 | PASS 8 ] ══ Failed tests ════════════════════════════════════════════════════════════════ ── Error ('test-glcm_textures.R:45:3'): glcm_textures ep works ───────────────── Error: [rast] empty srs Backtrace: ▆ 1. ├─terra::rast(...) at test-glcm_textures.R:45:3 2. └─terra::rast(...) 3. └─terra (local) .local(x, ...) 4. ├─terra::rast(nrows = nrow(x), ncols = ncol(x), crs = crs, extent = extent) 5. └─terra::rast(nrows = nrow(x), ncols = ncol(x), crs = crs, extent = extent) 6. └─terra (local) .local(x = x, ...) 7. └─terra:::new_rast(...) 8. └─terra:::messages(r, "rast") 9. └─terra:::error(f, x@pntr$getError()) ── Error ('test-glcm_textures.R:54:3'): glcm_textures er works ───────────────── Error: [rast] empty srs Backtrace: ▆ 1. ├─terra::rast(...) at test-glcm_textures.R:54:3 2. └─terra::rast(...) 3. └─terra (local) .local(x, ...) 4. ├─terra::rast(nrows = nrow(x), ncols = ncol(x), crs = crs, extent = extent) 5. └─terra::rast(nrows = nrow(x), ncols = ncol(x), crs = crs, extent = extent) 6. └─terra (local) .local(x = x, ...) 7. └─terra:::new_rast(...) 8. └─terra:::messages(r, "rast") 9. └─terra:::error(f, x@pntr$getError()) ── Error ('test-glcm_textures.R:63:3'): glcm_textures NA handling and ordering of metrics works ── Error: [rast] empty srs Backtrace: ▆ 1. ├─terra::rast(...) at test-glcm_textures.R:63:3 2. └─terra::rast(...) 3. └─terra (local) .local(x, ...) 4. ├─terra::rast(nrows = nrow(x), ncols = ncol(x), crs = crs, extent = extent) 5. └─terra::rast(nrows = nrow(x), ncols = ncol(x), crs = crs, extent = extent) 6. └─terra (local) .local(x = x, ...) 7. └─terra:::new_rast(...) 8. └─terra:::messages(r, "rast") 9. └─terra:::error(f, x@pntr$getError()) ── Error ('test-glcm_textures.R:75:3'): glcm_textures NA handling and ordering of metrics works ── Error: [rast] empty srs Backtrace: ▆ 1. ├─terra::rast(...) at test-glcm_textures.R:75:3 2. └─terra::rast(...) 3. └─terra (local) .local(x, ...) 4. ├─terra::rast(nrows = nrow(x), ncols = ncol(x), crs = crs, extent = extent) 5. └─terra::rast(nrows = nrow(x), ncols = ncol(x), crs = crs, extent = extent) 6. └─terra (local) .local(x = x, ...) 7. └─terra:::new_rast(...) 8. └─terra:::messages(r, "rast") 9. └─terra:::error(f, x@pntr$getError()) [ FAIL 4 | WARN 4 | SKIP 0 | PASS 8 ] Error: ! Test failures. Execution halted Flavor: r-oldrel-windows-x86_64

Version: 0.6.3
Check: re-building of vignette outputs
Result: ERROR Error(s) in re-building vignettes: --- re-building 'README.Rmd' using rmarkdown <!-- badges: start --> [](https://github.com/ailich/GLCMTextures/actions/workflows/R-CMD-check.yaml) [](https://cran.r-project.org/package=GLCMTextures) [](https://www.gnu.org/licenses/gpl-3.0) [](https://zenodo.org/badge/latestdoi/299630902) <!-- badges: end --> Please cite as Ilich, Alexander R. 2020. "GLCMTextures.", https://doi.org/10.5281/zenodo.4310186. https://github.com/ailich/GLCMTextures. # GLCMTextures ## Purpose This R package calculates the most common gray-level co-occurrence matrix (GLCM) texture metrics used for spatial analysis (Hall-Beyer 2017). It interfaces with C++ via the Rcpp and RcppArmadillo packages for increased speed. Texture metrics are calculated using a symmetric gray level co-occurence matrix (GLCM), meaning that each pixel is counted as a focal and neighboring pixel. For more details on how a symmetric GLCM is constructed, I highly recommend checking out Dr. Mryka Hall-Beyer's [GLCM texture tutorial](https://prism.ucalgary.ca/bitstream/handle/1880/51900/texture%20tutorial%20v%203_0%20180206.pdf). ## Motivation When comparing results across different software that calculate GLCM texture metrics, there are inconsistencies among results. This package is meant to provide a clearly documented implementation of GLCM texture metrics that gives the user control over key parameters to make it clear to the user exactly what they are calculating. As such, the formulas for each texture metric are provided, different shifts can be specified, the user can decide how to handle NA values, and the user gets control over how and if the data should be quantized. ## Install and Load Package The package can be installed from CRAN using `install.packages("GLCMTextures")` or the development version can be installed from github using the code `remotes::install_github("ailich/GLCMTextures")`. If you are using Windows, you may need to install Rtools using the instructions found [here](https://cran.r-project.org/bin/windows/Rtools/)). To install from github you must already have the remotes package installed, which can be installed using `install.packages("remotes")` This package relies on the `terra` package for handling of spatial raster data. ## Specifying the Relationship Between Focal and Neighbor Pixels The convention for specifying the direction of the neighboring pixel (the shift) is shown in the image below. The blue pixel in the center is treated as the focal pixel in the example. Shifts are specified as `c(x_shift, y_shift)`. So, a shift of `c(1,0)` refers to a the neighboring pixel being 1 to the right and 0 upwards of the focal pixel. Since a symmetric GLCM is created, this means each pixel is counted as both a focal and a neighboring pixel, so it also tabulates the shift in the opposite direction `c(-1,0)`, which is the dotted blue line. Therefore, these two shifts will produce equivalent results. Although neighboring pixels are typically considered as those one away in a given direction, the shift value can be specified as any integer value.  ## Available Metrics There are 8 metrics than can be calculated by this package that can be divided into 3 groups: the contrast group, the orderliness group, and the descriptive statistics group (Hall-Beyer 2017). The formulas provided below are from Hall-Beyer (2017). N = Number of rows or columns in the GLCM (Equal to the number of gray levels) i = row indices of the GLCM matrix (equal to gray level of reference cell) j = column indices of the GLCM matrix (equal to gray level of neighboring cell) P~i,j~ = Probability (relative frequency) of neighboring cells having gray levels i & j ### Contrast Group $$\text{GLCM Contrast} = \sum_{i,j=0}^{N-1} {P_{i,j}(i-j)^2}$$ $$\text{GLCM Dissimilarity} = \sum_{i,j=0}^{N-1} {P_{i,j}|i-j|}$$ $$\text{GLCM Homogeneity} = \sum_{i,j=0}^{N-1} \frac{P_{i,j}}{1+(i-j)^2}$$ ### Orderliness Group $$\text{GLCM Angular Second Moment (ASM)} = \sum_{i,j=0}^{N-1} {P_{i,j}^2}$$ $$\text{GLCM Entropy} = \sum_{i,j=0}^{N-1} {P_{i,j}[-ln(P_{i,j})]} \text{ where } 0*ln(0)=0$$ ### Descriptive Statistics Group $$\text{GLCM Mean} (\mu) = \sum_{i,j=0}^{N-1} i(P_{i,j})$$ $$\text{GLCM Variance} (\sigma^2) = \sum_{i,j=0}^{N-1} P_{i,j}(i-\mu)^2$$ $$\text{GLCM Correlation} = \sum_{i,j=0}^{N-1} {P_{i,j} \frac{(i-\mu)(j-\mu)}{\sigma^2}}$$ ## Tutorial Load packages ``` r library(GLCMTextures) #Load GLCMTextures package ``` See package help page ``` r help(package="GLCMTextures") ``` ### Test Matrix Before conducting texture calculations on entire raster data sets, we will work with a small matrix. ``` r test_matrix<- matrix(data=c(2,0,1,3,0,0,0,3,2), nrow = 3, ncol=3) print(test_matrix) #> [,1] [,2] [,3] #> [1,] 2 3 0 #> [2,] 0 0 3 #> [3,] 1 0 2 ``` This test matrix has 3 rows and 3 columns and contains values from 0-3 (4 gray levels). We can use the`make_glcm` function to create a normalized symmetric GLCM. A GLCM is a tabulation of counts and has the dimensions of the number of gray levels. The GLCM is initialized with all zeros and then we add as we tabulate counts. **Initialzed GLCM** ``` #> [,1] [,2] [,3] [,4] #> [1,] 0 0 0 0 #> [2,] 0 0 0 0 #> [3,] 0 0 0 0 #> [4,] 0 0 0 0 ``` The row and column number refers to the gray value of the focal and neighboring pixel (Since gray levels start at a value of 0, the row/column number is 1 larger than the corresponding gray level).We will use a shift of `c(1,0)` meaning that the neighboring pixels is the pixel directly to the right of the focal pixel. We start in the top left corner and we can see that we have a 2 as the focal value and a 3 as the neighboring value directly to the right, so we add 1 to the corresponding position in the GLCM which is row 3 (2+1)/column 4(3+1). Since we would like to create a symmetric GLCM where each pixel is treated as both a focal and neighbor value, we also add to row 4/column 3. ``` #> [,1] [,2] [,3] [,4] #> [1,] 0 0 0 0 #> [2,] 0 0 0 0 #> [3,] 0 0 0 1 #> [4,] 0 0 1 0 ``` We then continue this process throughout the whole image, moving right to the next focal pixel, and down to start the next row when a given row is completed. The resulting GLCM is a square matrix of counts that is symmetric about the diagonal. ``` r horizontal_glcm<- make_glcm(test_matrix, n_levels = 4, shift = c(1,0), normalize = FALSE) horizontal_glcm #> [,1] [,2] [,3] [,4] #> [1,] 2 1 1 2 #> [2,] 1 0 0 0 #> [3,] 1 0 0 1 #> [4,] 2 0 1 0 ``` Once we have finished tabulating all the counts we "normailize" the GLCM by dividing the each element by the sum of all the counts to get relative frequencies or probabilities that a given pixel of value i occurs next to a pixel of value j. The values in a normalized GLCM will therefore sum to 1. ``` r horizontal_glcm<- horizontal_glcm/sum(horizontal_glcm) horizontal_glcm #> [,1] [,2] [,3] [,4] #> [1,] 0.16666667 0.08333333 0.08333333 0.16666667 #> [2,] 0.08333333 0.00000000 0.00000000 0.00000000 #> [3,] 0.08333333 0.00000000 0.00000000 0.08333333 #> [4,] 0.16666667 0.00000000 0.08333333 0.00000000 ``` This could be accomplished in one line of code by setting the argument `normalize=TRUE` which is the default. ``` r make_glcm(test_matrix, n_levels = 4, shift = c(1,0), normalize = TRUE) #> [,1] [,2] [,3] [,4] #> [1,] 0.16666667 0.08333333 0.08333333 0.16666667 #> [2,] 0.08333333 0.00000000 0.00000000 0.00000000 #> [3,] 0.08333333 0.00000000 0.00000000 0.08333333 #> [4,] 0.16666667 0.00000000 0.08333333 0.00000000 ``` You may have noticed that pixels in the last column of the test matrix did not have a neighboring pixel to the right, so you would not tabulate any counts in those cases; however, this is precisely why we tabulate a symmetrical GLCM as these pixels do have neighbors to the left. Also, note that although the original matrix was 3x3, the GLCM is 4x4 because the size of the GLCM corresponds to the number of gray levels, not the size of the input matrix. Once the GLCM has been constructed, we can use this to calculate texture metrics using the`glcm_metrics` function to calculate the GLCM texture metrics. ``` r glcm_metrics(horizontal_glcm) #> glcm_contrast glcm_dissimilarity glcm_homogeneity glcm_ASM #> 4.000000 1.666667 0.400000 0.125000 #> glcm_entropy glcm_mean glcm_variance glcm_correlation #> 2.138333 1.166667 1.638889 -0.220339 ``` ### Raster Data Now we can move from calculating a single value of a given texture metric to calculating raster surfaces of texture metrics. Quitting from ./../man/fragments/README_Frag.Rmd:140-143 [elevation] ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ <error/rlang_error> Error: ! [rast] empty srs --- Backtrace: ▆ 1. ├─terra::rast(...) 2. └─terra::rast(...) 3. └─terra (local) .local(x, ...) 4. ├─terra::rast(nrows = nrow(x), ncols = ncol(x), crs = crs, extent = extent) 5. └─terra::rast(nrows = nrow(x), ncols = ncol(x), crs = crs, extent = extent) 6. └─terra (local) .local(x = x, ...) 7. └─terra:::new_rast(...) 8. └─terra:::messages(r, "rast") 9. └─terra:::error(f, x@pntr$getError()) ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Quitting from README.Rmd:23-24 [unnamed-chunk-2] ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ <error/rlang_error> Error: ! [rast] empty srs --- Backtrace: ▆ 1. ├─terra::rast(...) 2. └─terra::rast(...) 3. └─terra (local) .local(x, ...) 4. ├─terra::rast(nrows = nrow(x), ncols = ncol(x), crs = crs, extent = extent) 5. └─terra::rast(nrows = nrow(x), ncols = ncol(x), crs = crs, extent = extent) 6. └─terra (local) .local(x = x, ...) 7. └─terra:::new_rast(...) 8. └─terra:::messages(r, "rast") 9. └─terra:::error(f, x@pntr$getError()) ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Error: processing vignette 'README.Rmd' failed with diagnostics: [rast] empty srs --- failed re-building 'README.Rmd' SUMMARY: processing the following file failed: 'README.Rmd' Error: Vignette re-building failed. Execution halted Flavor: r-oldrel-windows-x86_64

Package MultiscaleDTM

Current CRAN status: ERROR: 2, NOTE: 2, OK: 9

Version: 1.0
Check: re-building of vignette outputs
Result: ERROR Error(s) in re-building vignettes: ... --- re-building ‘README.Rmd’ using rmarkdown # MultiscaleDTM <!-- badges: start --> [](https://github.com/ailich/MultiscaleDTM/actions/workflows/R-CMD-check.yaml) [](https://cran.r-project.org/package=MultiscaleDTM) [](https://www.gnu.org/licenses/gpl-3.0) [](https://zenodo.org/badge/latestdoi/353158828) <!-- badges: end --> Please cite as Ilich, A. R., Misiuk, B., Lecours, V., & Murawski, S. A. (2023). MultiscaleDTM: An open-source R package for multiscale geomorphometric analysis. _Transactions in GIS_, *27*(4). <https://doi.org/10.1111/tgis.13067> ## Purpose This R package calculates multi-scale geomorphometric terrain attributes from regularly gridded digital terrain models (DTM; i.e. elevation or bathymetry rasters) via a specified window as described in [Ilich et al. (2023)](http://doi.org/10.1111/tgis.13067). <img src="../man/figures/kwindow.png" width="50%"> Figure adapted from Wilson et al. (2007) ## Install and Load Package The package can be installed from CRAN using `install.packages("MultiscaleDTM")` or the development version can be installed from github using the code `remotes::install_github("ailich/MultiscaleDTM")`. If you are using Windows, you may need to install Rtools using the instructions found [here](https://cran.r-project.org/bin/windows/Rtools/)). To install from github you must already have the remotes package installed, which can be installed using `install.packages("remotes")` This package relies on the `terra` package for handling of spatial raster data. ## Main Functions ### Slope, Aspect and Curvature - `SlpAsp` - Calculates multi-scale slope and aspect using a modified version of the algorithm from Misiuk et al (2021) which extends classical formulations of slope restricted to a 3x3 window (Fleming and Hoffer, 1979; Horn et al., 1981; Ritter, 1987). This algorithm only considers a subset of cells within the focal window, specifically the four cells on the edge of the focal window directly up, down, left, and right of the focal cell for the "rook" method, an additional four corner cells for the "queen" method, or all edge cells for the "boundary" method. <img src="../man/figures/SlpAsp.png" width="70%"> - `DirSlp` - Calculates multi-scale slope in a specified direction. - `Qfit` - Calculates slope, aspect, curvature, and morphometric features by fitting a quadratic surface to the focal window using ordinary least squares using the equation shown below where a-f are regression parameters, Z is the elevation/depth, X is the east/west coordinates in the focal window relative to the focal cell, and Y is the north/south coordinates in the focal window relative to the focal cell (Evans, 1980; Wilson et al., 2007; Wood, 1996). The morphometric features algorithm has been modified to use more robust measures of curvature based on the suggestions of Minár et al. (2020). $$ Z = aX^2 + bY^2 +cXY+ dX +eY +f $$ <img src="../man/figures/Qfit_annotated.png" width="70%"> Figure adapted from Walbridge et al., (2018) - `Pfit` - Calculates multiscale slope and aspect aspect by fitting a planar surface to the focal window using ordinary least squares. This will provide equivalent results to using a quadratic fit (Jones, 1998) but is less computationally expensive. ### Roughness - `VRM` - Vector ruggedness measure (Sappington et al. 2007) quantifies roughness by measuring the dispersion of vectors normal to the terrain surface. This is accomplished by calculating the local (3 x 3 cell) slope and aspect, and constructing unit vectors normal to each cell in the DTM. These unit vectors are then decomposed into their corresponding x, y, and z components (i.e. the x, y, and z coordinates of the head of the vector relative to its origin) and used in the following equation (note: N is the number of cells in the window). VRM ranges from zero to one, representing completely smooth to rough surfaces, respectively. <img src="../man/figures/VRM_annotated.png" width="70%"> Figure adapted from Sappington et al. (2007) <img src="../man/figures/VRM_focal.png" width="70%"> Figure adapted from Habib (2021) $$ \text{VRM} = 1- \frac{\sqrt{\bigg(\sum x\bigg)^2+\bigg(\sum y\bigg)^2+\bigg(\sum z\bigg)^2}}{N} $$ $$ x = sin(\text{slope})*sin(\text{aspect}) $$ $$ y=sin(\text{slope})*cos(\text{aspect}) $$ $$ z=cos(\text{slope}) $$ - `SAPA` - Calculates the Surface Area to Planar Area (Jenness, 2004). Rougher surfaces will have a greater surface area to planar area ratio, and perfectly smooth surfaces will have a value of 1. This is a 3D analog to the classical "chain-and-tape" method, which calculates roughness as the ratio of the contoured distance (chain length) and linear distance (tape measure distance; Risk, 1972). Additionally, planar area can be corrected for slope by dividing the product of the x and y resolution by the cosine of slope (Du Preez 2015). The metric by Jenness (2004) and De Preez (2015) works at the per cell level (1x1 cell). This function generalizes this method to multiple scales by summing the surface areas within the focal window and adjusting the planar area of the focal window using multi-scale slope (Ilich et al., 2023). - `SurfaceArea` - Calculate the surface area of each grid cell (Jenness, 2004). This is accomplished by connecting a focal cell to its immediate neighbors to create 8 large triangles. These large triangles are then trimmed back to the extent of the focal cell using the principle of similar triangles, and then the area of those 8 smaller triangles are calculated and summed to estimate the surface area of the focal pixel. This is used within `SAPA`. <img src="../man/figures/chaintape.png" width="90%"> Figure adapted from Friedman et al. (2012) and created with BioRender.com. <img src="../man/figures/SAPA_annotated.png" width="70%"> Figure adapted from Jenness (2004) - `AdjSD`- This roughness metric modifies the standard deviation of elevation/bathymetry to account for slope (Ilich et al., 2023). It does this by first fitting a plane to the data in the focal window using ordinary least squares, and then extracting the residuals, and then calculating the standard deviation of the residuals within the focal window. <img src="../man/figures/adj_sd.png" width="80%"> - `RIE` - Calculates the Roughness Index-Elevation which quantifies the standard deviation of residual topography (Cavalli et al., 2008). This measure is conceptually similar to `AdjSD` but rather than fitting a plane and extracting residuals for the entire focal window, residual topography is calculated as the focal pixel minus the focal mean. Then the local standard deviation is calculated from this residual topography using a focal filter. <img src="../man/figures/RIE.png" width="80%"> Figure adapted from Cavalli et al. (2008) ### Relative Position Relative position represents whether an area is a local high or low in relation to a reference height. It is calculated as the value of the focal cell minus the value of a reference elevation which is often the mean of included values in the focal window but could also be other functions such as the median, min, or max of included values. Positive values indicate local topographic highs and negative values indicate lows.Relative Position can be expressed in units of the input DTM raster or can standardized relative to the local topography by dividing by the standard deviation or range of included elevation values in the focal window. - `RelPos` - A flexible and general purpose function to calculate relative position using a rectangular, circular, annulus, or custom shaped focal window and various functions of the included values as the reference height All other relative position functions are calls to `RelPos` with different default parameter values. - `TPI` - Topographic Position Index (Weiss, 2001) is the difference between the value of a focal cell and the mean of the surrounding cells (i.e. the central cell is excluded from focal opertaions) within a rectangular or circular focal window. - `DMV` - Difference from Mean Value (Lecours et al., 2017; Wilson, and Gallant, 2000) is the difference between the value of a focal cell and the mean of all cells (i.e. including the focal cell) in a rectangular or circular focal window. - `BPI` - Bathymetric Position Index (Lundblad et al., 2006) is the difference between the value of a focal cell and the mean of the surrounding cells contained within an annulus shaped window. Since an annulus shaped window is used, it requires an inner and outer radius to be specified. Although the name contains "bathymetric," that is due to the context in which it was proposed, and is equally applicable to terrestrial elevation data. <img src="../man/figures/WindowShapes.png" width="100%"> Examples of different focal window shapes. Shown are a 13 x 13 cell rectangular window (left), a circular window with a radius of six cells (center), and an annulus window with an inner radius of four cells and an outer radius of six cells (right). ## Tutorial In this tutorial we will calculate most terrain attributes using a 5 x 5 cell rectangular window; however, any rectangular odd numbered window size could be used. Window size is specified using the `w` parameter. Rectangular window sizes are specified with a vector of length 2 as `c(n_rows, n_cols)`. If a single number is provided it will be used for both the number of rows and columns. Functions that calculate relative position currently support other focal window shapes. In those examples, we will additionally calculate the measures using a circular focal window with a radius of 2 cells, an annulus window with an inner radius of 4 and an outer radius of 6 cells, and a custom focal window. Circular windows are specified by a single number representing the radius, annulus windows are specified with a vector of length 2 of `c(inner_radius, outer_radius)`, and custom windows are specified by a matrix with values showing which data to include (1's) and which data to exclude (NA's). **Load packages** ``` r library(MultiscaleDTM) #Load MultiscaleDTM package ``` **See package help page** ``` r help(package="MultiscaleDTM") ``` **Read in Data** Typically raster data would be loaded into R via `terra` package's `rast` function. We will use instead load a georeferenced version of R's built in volcano data set using the `erupt` function. ``` r r<- erupt() ``` Quitting from ./../man/fragments/README_Frag.Rmd:141-146 [Topo] ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ <error/rlang_error> Error in `library()`: ! there is no package called 'tmap' --- Backtrace: ▆ 1. └─base::library(tmap) ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Quitting from README.Rmd:22-23 [unnamed-chunk-2] ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ <error/rlang_error> Error in `library()`: ! there is no package called 'tmap' --- Backtrace: ▆ 1. └─base::library(tmap) ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Error: processing vignette 'README.Rmd' failed with diagnostics: there is no package called 'tmap' --- failed re-building ‘README.Rmd’ SUMMARY: processing the following file failed: ‘README.Rmd’ Error: Vignette re-building failed. Execution halted Flavor: r-patched-linux-x86_64

Version: 1.0.1
Check: installed package size
Result: NOTE installed size is 8.3Mb sub-directories of 1Mb or more: doc 2.7Mb help 2.3Mb libs 2.2Mb Flavors: r-oldrel-macos-arm64, r-oldrel-macos-x86_64, r-oldrel-windows-x86_64

Version: 1.0.1
Check: examples
Result: ERROR Running examples in 'MultiscaleDTM-Ex.R' failed The error most likely occurred in: > ### Name: AdjSD > ### Title: Calculates standard deviation of bathymetry (a measure of > ### rugosity) adjusted for slope > ### Aliases: AdjSD > > ### ** Examples > > r<- erupt() Warning: PROJ: proj_create_from_database: Cannot find proj.db (GDAL error 1) Error: [rast] empty srs Execution halted Flavor: r-oldrel-windows-x86_64

Version: 1.0.1
Check: tests
Result: ERROR Running 'testthat.R' [11s] Running the tests in 'tests/testthat.R' failed. Complete output: > # This file is part of the standard setup for testthat. > # It is recommended that you do not modify it. > # > # Where should you do additional test configuration? > # Learn more about the roles of various files in: > # * https://r-pkgs.org/testing-design.html#sec-tests-files-overview > # * https://testthat.r-lib.org/articles/special-files.html > > library(testthat) > library(MultiscaleDTM) Loading required package: terra terra 1.8.86 Attaching package: 'terra' The following objects are masked from 'package:testthat': compare, describe > > test_check("MultiscaleDTM") Saving _problems/test-functions-1.R [ FAIL 1 | WARN 1 | SKIP 0 | PASS 0 ] ══ Failed tests ════════════════════════════════════════════════════════════════ ── Error ('test-functions.R:1:1'): (code run outside of `test_that()`) ───────── Error: [rast] empty srs Backtrace: ▆ 1. └─MultiscaleDTM::erupt() at test-functions.R:1:1 2. ├─terra::rast(...) 3. └─terra::rast(...) 4. └─terra (local) .local(x, ...) 5. ├─terra::rast(nrows = nrow(x), ncols = ncol(x), crs = crs, extent = extent) 6. └─terra::rast(nrows = nrow(x), ncols = ncol(x), crs = crs, extent = extent) 7. └─terra (local) .local(x = x, ...) 8. └─terra:::new_rast(...) 9. └─terra:::messages(r, "rast") 10. └─terra:::error(f, x@pntr$getError()) [ FAIL 1 | WARN 1 | SKIP 0 | PASS 0 ] Error: ! Test failures. Execution halted Flavor: r-oldrel-windows-x86_64

Version: 1.0.1
Check: re-building of vignette outputs
Result: ERROR Error(s) in re-building vignettes: --- re-building 'README.Rmd' using rmarkdown # MultiscaleDTM <!-- badges: start --> [](https://github.com/ailich/MultiscaleDTM/actions/workflows/R-CMD-check.yaml) [](https://cran.r-project.org/package=MultiscaleDTM) [](https://www.gnu.org/licenses/gpl-3.0) [](https://zenodo.org/badge/latestdoi/353158828) <!-- badges: end --> Please cite as Ilich, A. R., Misiuk, B., Lecours, V., & Murawski, S. A. (2023). MultiscaleDTM: An open-source R package for multiscale geomorphometric analysis. _Transactions in GIS_, *27*(4). <https://doi.org/10.1111/tgis.13067> ## Purpose This R package calculates multi-scale geomorphometric terrain attributes from regularly gridded digital terrain models (DTM; i.e. elevation or bathymetry rasters) via a specified window as described in [Ilich et al. (2023)](http://doi.org/10.1111/tgis.13067). <img src="../man/figures/kwindow.png" width="50%"> Figure adapted from Wilson et al. (2007) ## Install and Load Package The package can be installed from CRAN using `install.packages("MultiscaleDTM")` or the development version can be installed from github using the code `remotes::install_github("ailich/MultiscaleDTM")`. If you are using Windows, you may need to install Rtools using the instructions found [here](https://cran.r-project.org/bin/windows/Rtools/)). To install from github you must already have the remotes package installed, which can be installed using `install.packages("remotes")` This package relies on the `terra` package for handling of spatial raster data. ## Main Functions ### Slope, Aspect and Curvature - `SlpAsp` - Calculates multi-scale slope and aspect using a modified version of the algorithm from Misiuk et al (2021) which extends classical formulations of slope restricted to a 3x3 window (Fleming and Hoffer, 1979; Horn et al., 1981; Ritter, 1987). This algorithm only considers a subset of cells within the focal window, specifically the four cells on the edge of the focal window directly up, down, left, and right of the focal cell for the "rook" method, an additional four corner cells for the "queen" method, or all edge cells for the "boundary" method. <img src="../man/figures/SlpAsp.png" width="70%"> - `DirSlp` - Calculates multi-scale slope in a specified direction. - `Qfit` - Calculates slope, aspect, curvature, and morphometric features by fitting a quadratic surface to the focal window using ordinary least squares using the equation shown below where a-f are regression parameters, Z is the elevation/depth, X is the east/west coordinates in the focal window relative to the focal cell, and Y is the north/south coordinates in the focal window relative to the focal cell (Evans, 1980; Wilson et al., 2007; Wood, 1996). The morphometric features algorithm has been modified to use more robust measures of curvature based on the suggestions of Minár et al. (2020). $$ Z = aX^2 + bY^2 +cXY+ dX +eY +f $$ <img src="../man/figures/Qfit_annotated.png" width="70%"> Figure adapted from Walbridge et al., (2018) - `Pfit` - Calculates multiscale slope and aspect aspect by fitting a planar surface to the focal window using ordinary least squares. This will provide equivalent results to using a quadratic fit (Jones, 1998) but is less computationally expensive. ### Roughness - `VRM` - Vector ruggedness measure (Sappington et al. 2007) quantifies roughness by measuring the dispersion of vectors normal to the terrain surface. This is accomplished by calculating the local (3 x 3 cell) slope and aspect, and constructing unit vectors normal to each cell in the DTM. These unit vectors are then decomposed into their corresponding x, y, and z components (i.e. the x, y, and z coordinates of the head of the vector relative to its origin) and used in the following equation (note: N is the number of cells in the window). VRM ranges from zero to one, representing completely smooth to rough surfaces, respectively. <img src="../man/figures/VRM_annotated.png" width="70%"> Figure adapted from Sappington et al. (2007) <img src="../man/figures/VRM_focal.png" width="70%"> Figure adapted from Habib (2021) $$ \text{VRM} = 1- \frac{\sqrt{\bigg(\sum x\bigg)^2+\bigg(\sum y\bigg)^2+\bigg(\sum z\bigg)^2}}{N} $$ $$ x = sin(\text{slope})*sin(\text{aspect}) $$ $$ y=sin(\text{slope})*cos(\text{aspect}) $$ $$ z=cos(\text{slope}) $$ - `SAPA` - Calculates the Surface Area to Planar Area (Jenness, 2004). Rougher surfaces will have a greater surface area to planar area ratio, and perfectly smooth surfaces will have a value of 1. This is a 3D analog to the classical "chain-and-tape" method, which calculates roughness as the ratio of the contoured distance (chain length) and linear distance (tape measure distance; Risk, 1972). Additionally, planar area can be corrected for slope by dividing the product of the x and y resolution by the cosine of slope (Du Preez 2015). The metric by Jenness (2004) and De Preez (2015) works at the per cell level (1x1 cell). This function generalizes this method to multiple scales by summing the surface areas within the focal window and adjusting the planar area of the focal window using multi-scale slope (Ilich et al., 2023). - `SurfaceArea` - Calculate the surface area of each grid cell (Jenness, 2004). This is accomplished by connecting a focal cell to its immediate neighbors to create 8 large triangles. These large triangles are then trimmed back to the extent of the focal cell using the principle of similar triangles, and then the area of those 8 smaller triangles are calculated and summed to estimate the surface area of the focal pixel. This is used within `SAPA`. <img src="../man/figures/chaintape.png" width="90%"> Figure adapted from Friedman et al. (2012) and created with BioRender.com. <img src="../man/figures/SAPA_annotated.png" width="70%"> Figure adapted from Jenness (2004) - `AdjSD`- This roughness metric modifies the standard deviation of elevation/bathymetry to account for slope (Ilich et al., 2023). It does this by first fitting a plane to the data in the focal window using ordinary least squares, and then extracting the residuals, and then calculating the standard deviation of the residuals within the focal window. <img src="../man/figures/adj_sd.png" width="80%"> - `RIE` - Calculates the Roughness Index-Elevation which quantifies the standard deviation of residual topography (Cavalli et al., 2008). This measure is conceptually similar to `AdjSD` but rather than fitting a plane and extracting residuals for the entire focal window, residual topography is calculated as the focal pixel minus the focal mean. Then the local standard deviation is calculated from this residual topography using a focal filter. <img src="../man/figures/RIE.png" width="80%"> Figure adapted from Cavalli et al. (2008) ### Relative Position Relative position represents whether an area is a local high or low in relation to a reference height. It is calculated as the value of the focal cell minus the value of a reference elevation which is often the mean of included values in the focal window but could also be other functions such as the median, min, or max of included values. Positive values indicate local topographic highs and negative values indicate lows.Relative Position can be expressed in units of the input DTM raster or can standardized relative to the local topography by dividing by the standard deviation or range of included elevation values in the focal window. - `RelPos` - A flexible and general purpose function to calculate relative position using a rectangular, circular, annulus, or custom shaped focal window and various functions of the included values as the reference height All other relative position functions are calls to `RelPos` with different default parameter values. - `TPI` - Topographic Position Index (Weiss, 2001) is the difference between the value of a focal cell and the mean of the surrounding cells (i.e. the central cell is excluded from focal opertaions) within a rectangular or circular focal window. - `DMV` - Difference from Mean Value (Lecours et al., 2017; Wilson, and Gallant, 2000) is the difference between the value of a focal cell and the mean of all cells (i.e. including the focal cell) in a rectangular or circular focal window. - `BPI` - Bathymetric Position Index (Lundblad et al., 2006) is the difference between the value of a focal cell and the mean of the surrounding cells contained within an annulus shaped window. Since an annulus shaped window is used, it requires an inner and outer radius to be specified. Although the name contains "bathymetric," that is due to the context in which it was proposed, and is equally applicable to terrestrial elevation data. <img src="../man/figures/WindowShapes.png" width="100%"> Examples of different focal window shapes. Shown are a 13 x 13 cell rectangular window (left), a circular window with a radius of six cells (center), and an annulus window with an inner radius of four cells and an outer radius of six cells (right). ## Tutorial In this tutorial we will calculate most terrain attributes using a 5 x 5 cell rectangular window; however, any rectangular odd numbered window size could be used. Window size is specified using the `w` parameter. Rectangular window sizes are specified with a vector of length 2 as `c(n_rows, n_cols)`. If a single number is provided it will be used for both the number of rows and columns. Functions that calculate relative position currently support other focal window shapes. In those examples, we will additionally calculate the measures using a circular focal window with a radius of 2 cells, an annulus window with an inner radius of 4 and an outer radius of 6 cells, and a custom focal window. Circular windows are specified by a single number representing the radius, annulus windows are specified with a vector of length 2 of `c(inner_radius, outer_radius)`, and custom windows are specified by a matrix with values showing which data to include (1's) and which data to exclude (NA's). **Load packages** ``` r library(MultiscaleDTM) #Load MultiscaleDTM package ``` **See package help page** ``` r help(package="MultiscaleDTM") ``` **Read in Data** Typically raster data would be loaded into R via `terra` package's `rast` function. We will use instead load a georeferenced version of R's built in volcano data set using the `erupt` function. Quitting from ./../man/fragments/README_Frag.Rmd:150-152 [unnamed-chunk-8] ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ <error/rlang_error> Error: ! [rast] empty srs --- Backtrace: ▆ 1. └─MultiscaleDTM::erupt() 2. ├─terra::rast(...) 3. └─terra::rast(...) 4. └─terra (local) .local(x, ...) 5. ├─terra::rast(nrows = nrow(x), ncols = ncol(x), crs = crs, extent = extent) 6. └─terra::rast(nrows = nrow(x), ncols = ncol(x), crs = crs, extent = extent) 7. └─terra (local) .local(x = x, ...) 8. └─terra:::new_rast(...) 9. └─terra:::messages(r, "rast") 10. └─terra:::error(f, x@pntr$getError()) ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Quitting from README.Rmd:22-23 [unnamed-chunk-2] ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ <error/rlang_error> Error: ! [rast] empty srs --- Backtrace: ▆ 1. └─MultiscaleDTM::erupt() 2. ├─terra::rast(...) 3. └─terra::rast(...) 4. └─terra (local) .local(x, ...) 5. ├─terra::rast(nrows = nrow(x), ncols = ncol(x), crs = crs, extent = extent) 6. └─terra::rast(nrows = nrow(x), ncols = ncol(x), crs = crs, extent = extent) 7. └─terra (local) .local(x = x, ...) 8. └─terra:::new_rast(...) 9. └─terra:::messages(r, "rast") 10. └─terra:::error(f, x@pntr$getError()) ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Error: processing vignette 'README.Rmd' failed with diagnostics: [rast] empty srs --- failed re-building 'README.Rmd' SUMMARY: processing the following file failed: 'README.Rmd' Error: Vignette re-building failed. Execution halted Flavor: r-oldrel-windows-x86_64