Generate neutral landscapes that match exactly the extent (non-NA cells) of another raster about landscapetools HOT 10 CLOSED

Hi André,

thanks for posting here! In your example it actually makes sense, because the masking "shrinks" the raster.

I somewhat faked the dataset you were sending via mail and attach an example that I hope better illustrates your problem:

library(NLMR)
library(raster)
#> Loading required package: sp

test_raster <- raster('test_raster.tif')

plot(test_raster, col=topo.colors(100))

# check how many cells are not NA, but what we actually need would be the 
# side length of the non-na cells - so your acctual data
sum(!is.na(test_raster[]))
#> [1] 169846

# Following from here is basically your approach André, this is what 
# mask does in the background (boiled down)
random_raster <- nlm_distancegradient(ncol = 588, nrow = 425, resolution = 1, origin = c(1,0,0,425))  # extent of test_raster, problem being that we also have values for the na-cells
plot(random_raster)

range(random_raster[], na.rm = TRUE)
#> [1] 0 1

# stripped down mask function
x <- getValues(random_raster)
mask <- getValues(test_raster)
x[is.na(mask) & !is.na(x)] <- NA
x <- setValues(random_raster, x)

# and we are missing values, since we cut pieces of
plot(x, col=topo.colors(100))

range(x[], na.rm = TRUE)
#> [1] 0.000000 0.976313

We obviously loose some values because amount of data in your input file is lower than the neutral landscape model we create. I think to solve that, we would have to reproject the input data in a euclidian projection (starting at 0,0) and simulate a NLM (same number of columns and rows as the trimmed input data) with the same properties. @Nowosad this is what I was talking about on gitter.

Here is the test_raster file:
test_raster.tif.tar.gz

^{Created on 2018-12-13 by the reprex package (v0.2.1)}

from landscapetools.

I cannot think of a good solution now. Maybe @mdsumner can think of something..

from landscapetools.

I appreciate your efforts. I still have not figure out a clean way to do this.

from landscapetools.

I don't fully understand yet, it seems to me that you get nlm values where you have data values by using mask. I can see in @marcosci's example that you probably want the raster aligned orthogonally rather than on an angle? But still I don't understand why this is a problem, could you expand on what happens next for why this is an issue? (That might help me understand).

BTW I also use this 0,ncol 0,nrow extent very often, though I never used "origin" or even noticed it before. I use setExent(0, ncol(x), 0, nrow(x)) - though that needs to be row/col flipped if input is a matrix.

It works brilliantly for dealing with tricky grids from ocean model output and I now use a "coords" 2-layer raster to keep track of the geography, and plot and even regrid the data - some of those tricks might help here.

from landscapetools.

@mdsumner my problem is actually a minor one. Still, when I create a neutral landscape model (nlm) based on the extent of a specific environmental raster layer, such as mine (if you want, you can use the code that I provide in the original post), the generated nlm layer will be larger in area than my the environmental raster that I used to set the extent. This distance nlm layer has a range of values between 0 and 1.

However, when I apply mask to this nlm layer to match exactly a given environmental raster layer , the values in the edges of the nlm layer are lost because the environmental raster that I used to mask is smaller in size. Therefore, the "masked" nlm layer will have values between something like 0.05 and 0.94 (not the full original range of 0-1), and I would like to preserve the full range of 0-1 values. The only solution that I was able to come up with was to make a second scaling to the masked nlm layer:

dist.grad2<-(dist.grad1-cellStats(dist.grad1,"min"))/(cellStats(dist.grad1,"max")-cellStats(dist.grad1,"min"))

I hope this explanation helps. Thank you very much for your help.

from landscapetools.

I think the biggest abstraction of the problem is that NLM are always orthogonal rectangles.
However, their purpose of being a spatial null model makes it necessary to simulate landscapes in the mask of real data, so the area where the extent of a real raster dataset (as here #15 (comment)) has data and is not NA.

To do that, we would need to find out the number of rows and columns of this dataset. In the example it's something like 420 x 530, but this also includes columns and rows with NA.
raster already has trim, but that only works if an entire row or column has only NA.

Is it possible to rectify/realing (I am obviously missing the correct word) projected data, to find out the extent of the area with values @mdsumner?

from landscapetools.

Oh I see, the range is lost. First cheat I think of is to rescale to 0, 1 after masking but does that also lose overall properties of the field?
I'll keep thinking about it

from landscapetools.

That loses the properties, unfortunately, which is why would need to realign to original dataset ... Thanks for thinking about it!

from landscapetools.

I still don't understand, you seem to be setting part of the top row to 1 ,i.e

## line 69 in nlm_distancegradient.R
 distancegradient[origin[1]:origin[2], origin[3]:origin[4]] <- 1

which is the top row, but counting across to nrows - that's why there's a white line partially at the top left of random_raster. It's definitely possible to isolate that rotated hull (always convex?) efficiently and determine distances, but I can't see what's supposed to actually happen here - we effectively get distances to that top row of cells, from 1:425 in this example. I don't understand why it's not just a rescale to 0,1 for the remaining values, but I appreciate I don't understand overall either.

So my shot in the dark is (and apology if I'm just off base):

r <- raster("C:\\Users\\michae_sum\\Downloads\\test_raster.tif")
ok <- !is.na(values(r))

## cells count from 1:ncell (top-left to right, then down)
xy <- xyFromCell(r, 1:425)  ## this is what `distancegradient[origin[1]:origin[2], origin[3]:origin[4]] <- 1` isolates 

dp <- distanceFromPoints(r, xy)

nlm <- setValues(r, NA_real_)
nlm[ok] <- scales::rescale(dp[ok], c(0, 1))
plot(nlm)

from landscapetools.

And for @poroc300's example

library(raster)
library(NLMR)

#create a raster
r<-raster(ncol=100, nrow=100, xmn=0,xmx=100,ymn=0,ymx=100)
values(r) <- runif(ncell(r), min=1,max=100)

#insert NA values in the edges
r[c(1:500, 9500:10000)]<-NA


ok <- !is.na(values(r))
xy <- xyFromCell(r, 1)  ## origin(0, 1, 0, 1) is the top left cell as best I can see
dp <- distanceFromPoints(r, xy)
nlm <- setValues(r, NA_real_)
nlm[ok] <- scales::rescale(dp[ok], c(0, 1))
plot(nlm)

Again, apologies if I miss the brief - just ignore if so ;)

from landscapetools.