Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Global Moran's I in autocor function does not correctly count neighbours #1592

Closed
NicholasBerryman opened this issue Aug 28, 2024 · 1 comment

Comments

@NicholasBerryman
Copy link

NicholasBerryman commented Aug 28, 2024

This line in the global Moran's I calculation (autocor.R line 156) will return NA for any cell with an adjacent NA cell, and thus does not correctly count neighbours for any cell with an adjacent NA.

W <- focal( zz, w=w, fun="sum")

In extreme cases this leads to a failure to calculate Global Moran's I (autocorr returns NA).
In less extreme cases this leads to slightly inflated values of I being returned.

Here's some code to replicate this issue in an extreme case:

exampleRast_base <- ifel((0==rast(diag(1,5,5))),NA,1) # Generate a diagonal raster -> an extreme case
exampleRast      <- mask(focal(exampleRast_base,3,'sum',na.rm = TRUE),exampleRast_base) # Adjust the raster to add some variance
x <- exampleRast
w <- matrix(c(1,1,1,1,0,1,1,1,1),3)
autocor(x, w, method="moran", global=TRUE) # -> Returns NA

And here's some code that I think will properly count neighbours if we replace line 156 with it. Notably, we need to apply a mask using the original raster after adding na.rm = TRUE, otherwise the autocorr function will produce slightly deflated values for I.

W <- mask(focal( zz, w=w, fun="sum", na.rm = TRUE),x)

I'm unsure whether a similar issue exists for the local or Geary calculations - I'm not particularly familiar with how they are calculated.

OS is Windows 10

> packageVersion("terra")
[1] ‘1.7.55’
> terra::gdal(lib="all")
   gdal    proj    geos 
"3.5.2" "8.2.1" "3.9.3" 
@rhijmans
Copy link
Member

Thank you very much and sorry for the long response time. Fixed as suggested.

netbsd-srcmastr pushed a commit to NetBSD/pkgsrc that referenced this issue Feb 8, 2025
# Version 1.8-15

## Bug Fixes

- `Readrds` Failed For Rasters With Timestep="Seconds"
  [#1711](Https://github.com/Rspatial/terra/issues/1711) by Pascal
  Oettli

- `divide<SpatVector>` always returned NULL
  [#1724](rspatial/terra#1724) by Márcia
  Barbosa

- `erase` failed in some cases
  [#1710](rspatial/terra#1710) by
  erkent-carb

## enhancements

- `bestMatch` now has argument "fun" to allow the use of different
  distance measures, and a <matrix> method

- `wrap` (and `writeRDS`) now captures varnames/longnames
  [#1719](rspatial/terra#1719) by Andrew
  Gene Brown

- improved raster metadata writing
  [#1714](rspatial/terra#1714) by Andrew Gene
  Brown

- `vect` and `writeVector` now properly read and write date and
  datetime
  data. [#1718](rspatial/terra#1718) by
  Andrew Gene Brown

- improved estimate of available memory on linux systems
  [#1506](rspatial/terra#1506) by Cedric
  Rossi

# version 1.8-10

Released 2025-01-13

## bug fixes

- `expanse<SpatRaster>(transform=TRUE)` crashed R when the crs was
  "local". [#1671](rspatial/terra#1671) by
  Michael Chirico

- `patches(values=TRUE)` wrapped around the edges
  [#1675](rspatial/terra#1675) by Michael
  Chirico

- `spin` now correctly handles spherical coordinates
  [#1576](rspatial/terra#1576) by jeanlobry

- `mosaic` sometimes crashed R
  [#1524](rspatial/terra#1524) by John
  Baums, Dave Klinges, and Hugh Graham.

- `spatSample` ignored argument "exp" when taking a random sample with
  na.rm=TRUE on a large raster
  [#1437](rspatial/terra#1437) by Babak
  Naimi

- `split<SpatVector,SpatVector>` did not work properly
  [#1619](rspatial/terra#1619) by Michael
  Sumner

- `autocor` improved handling of NA cells for global Moran computation
  [#1992](rspatial/terra#1592) by Nicholas
  Berryman

- `shade` is more
  memory-safe. [#1452](rspatial/terra#1452)
  by Francis van Oordt and Chris English

- fixed bug in `rasterize` revealed when using `crop(mask=TRUE)`
  [#1686](rspatial/terra#1686) by edixon1

- fixed `to_id = NA` bug in `nearest`
  [#1471](rspatial/terra#1471) by Mats
  Blomqvist

- better handling of date/unit
  [#1684](rspatial/terra#1684) and
  [#1688](rspatial/terra#1688) by Andrew
  Gene Brown

- `spatSample(method="regular")` on a raster with one column returned
  too many samples
  [#1362](rspatial/terra#1362) by Daniel R
  Schlaepfer


## enhancements

- `plot<SpatVector>` now uses the same default viridis color palette
  as `plot<SpatRaster>`
  [#1670](rspatial/terra#1670) by Márcia
  Barbosa

- `relate` now accepts relation="equals"
  [#1672](rspatial/terra#1672) by Krzysztof
  Dyba

- `init` now accepts additional arguments for function "fun"

- better handling of the 32 connections limitation set by the HDF4
  library [#1481](rspatial/terra#1481) by
  Dimitri Falk

- When using RStudio a once per session warning is given when using
  draw, sel or click
  [#1063](rspatial/terra#1063) by Sergei
  Kharchenko

- `distance<SpatRaster>` from lon and lat lines/polygons computes
  distance to the edges instead of the nodes
  [#1462](rspatial/terra#1462) by Derek
  Friend

- `distance<SpatVector,SpatVector>` now works for lon/lat data
  [#1615](rspatial/terra#1615) by Wencheng
  Lau-Medrano

- using overviews for faster plotting of COGs over http
  [#1353](rspatial/terra#1353) by Michael
  Sumner and [#1412](rspatial/terra#1412);
  and argument `plot(x, overview=)` to change the default behavior.

- `extract` with points is now faster for rasters accessed over http
  [#1504](rspatial/terra#1504) by Krzysztof
  Dyba

- `extract` with many points on very large rasters was slower in
  compared to doing the same with "raster" (which uses terra for
  that!) [#1584](rspatial/terra#1584) by
  Hassan Masoomi

- `merge` now has three alternative algorithms
  [1366](rspatial/terra#1366) by Hassan
  Masoomi and [#1650](rspatial/terra#1650)
  by Agustin Lobo


## new

- `$<SpatRaster>` can now be used to get a categorical SpatRaster with
  a different active category

- `scale_linear<SpatRaster>` method for linear scaling of cell values
  between a minimum and maximum value such as 0 and 1

- `distance` and related methods get argument "method" to choose the
  distance algorithm for lon/lat data
  [#1677](rspatial/terra#1677) by Márcia
  Barbosa

- `divide<SpatRaster>` and `divide<SpatVector>` methods

- `nseg` counts the number of segments in a SpatVector
  [#1647](rspatial/terra#1674) by Michael
  Chirico

- `extract` argument "search_radius" to extract values from the
  nearest raster cell that is not `NA`
  [#873](rspatial/terra#873) by
  matthewseanmarcus

- `combineLevels` to combine the levels of all layers
  [link](https://stackoverflow.com/questions/79340152/how-to-set-factor-levels-in-a-rast-stack-using-terra-when-different-levels-exi)
  on SO by Sam
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

3 participants
@NicholasBerryman @rhijmans and others