TWest
TWest

Reputation: 824

Change the value of a given number of raster cells based on another raster in R

I am working on a land-use/cover change simulation in R. I have two raster maps, one with 'land-use/cover classes' and another with 'deforestation risk info.' I am using the risk raster to identify the forest pixels (one of the land-use/cover classes) more likely to be deforested. So far I have a R code that works, here is an example that can be reproduced:

#create land-use/cover raster
real <- raster(nrows=25, ncols=25, vals=round(rnorm(625, 3), 0))
real[ real > 3 ] <- 3 #just so we end with 3 land-use/cover classes
real[ real < 1 ] <- 1 #just so we end with 3 land-use/cover classes
plot(real)

#function to create the deforestation risk raster created by someone else
rtnorm <- function(n, mean = 0, sd = 1, min = 0, max = 1) { 
  bounds <- pnorm(c(min, max), mean, sd)
  u <- runif(n, bounds[1], bounds[2])
  qnorm(u, mean, sd)
}

risk <- raster(nrows=25, ncols=25, vals=rtnorm(n = 625, .1, .9)) #deforestation risk raster
plot(risk)

#The actual analysis starts here:
forest <- real #read the land-use/cover raster
forest[ forest != 3 ] <- NA #all left is class 3 (let's say, forest) in the raster
plot(forest, col="green")

deforestation <- sum(forest, risk) #identify the forest pixels with highest risk
plot(deforestation)

deforestation[ deforestation <= 3.5 ] <- 0 #rule to deforest the forest pixels
deforestation[ deforestation > 0 ] <- 100 #mark the pixels that will be deforested
plot(deforestation)

simulation <- sum(real, deforestation)
simulation[ simulation > 100 ] <- 2 #I use this to mark the forest pixels to a different land-use/cover class
plot(simulation)

I would like to change the rule I am using to select the forest pixels that will be deforested (i.e., deforestation[ deforestation <= 3.5 ] <- 0 ). Instead of selecting a threshold value like 3.5, I wonder if I could set a specific number of forest pixels (e.g., 50) to be deforested, and then select the 50 forest pixels with the highest deforestation risk. I am completely clueless about how to do something like this in R, so any suggestions will be highly appreciated. Thank you.

Upvotes: 0

Views: 83

Answers (1)

Robert Hijmans
Robert Hijmans

Reputation: 47156

If the raster is not too large, you can get the 50 cells with the highest values as follows:

i <- order(values(deforestation), decreasing=TRUE)[1:50]

You can then use these like this

deforestation[] <- 0
deforestation[i] <- 100

Upvotes: 1

Related Questions