3

I have a raster:

r <- raster(ncol=10, nrow=10)
set.seed(0)
values(r) <- runif(ncell(r))

img1

From the raster I select the top 10% and change to binary:

r_10<-r[[1]]>=quantile(r,.90)

img2

From this subset raster r_10 all green pixels have the same value of 1. I would like to change these values, by identifying pixels or groups of pixels as objects and labeling every new object with a new ID. The new raster should have values like this example image: img3

Some objects can have multiple pixels, and they all should have the same object ID (like number 8).

How can I code this up in R? I thought to use some sort of edge detection, or Sobel filter, but cant figure it out.

Here is a similar post, not the same, but its in python, and I need to implement this in R.

Any alternative solutions are welcome.

kjtheron
  • 191
  • 1
  • 11
  • what does it mean by "objects"? All these green pixels have value of 1. When you label them, does the order matter? For example, can the number 1 be assigned to the bottom pixel instead of the top-left one? – uPhone Mar 05 '21 at 16:19
  • I use the term "object" very loosely. Here, I just want to reclassify groups of neighboring pixels to different values. The order of labeling them does not matter. As long as very "object" has its own unique ID. – kjtheron Mar 05 '21 at 16:24

1 Answers1

2

I am sure there are multiple ways to answer this questions (computer vision and GIS). Here is an GIS solution (found here) to the problem at hand:

# Create raster data
r <- raster(ncol=10, nrow=10)
set.seed(0)
values(r) <- runif(ncell(r))

# Select top 10% of highest values and convert to binary
r_10<-r[[1]]>=quantile(r,.90)
r_10[r_10==0]<-NA

# Vectorize
Vector_r_10<-rasterToPolygons(r_10)
plot(Vector_r_10)

# Add new Obj_ID class
Vector_r_10$Obj_ID<-1:nrow(Vector_r_10)

# Identify neighboring pixels
nb<-poly2nb(Vector_r_10)

# Create regions
create_regions <- function(data) {
  group <- rep(NA, length(data))
  group_val <- 0
  while(NA %in% group) {
    index <- min(which(is.na(group)))
    nb <- unlist(data[index])
    nb_value <- group[nb]
    is_na <- is.na(nb_value)
    if(sum(!is_na) != 0){
      prev_group <- nb_value[!is_na][1]
      group[index] <- prev_group
      group[nb[is_na]] <- prev_group
    } else {
      group_val <- group_val + 1
      group[index] <- group_val
      group[nb] <- group_val
    }
  }
  group
}
region<-create_regions(nb)

# Union on new regions
pol_rgn<-spCbind(Vector_r_10,region)
New_Vector_r_10<-unionSpatialPolygons(pol_rgn,region)
New_Vector_r_10<-as(New_Vector_r_10,"SpatialPolygonsDataFrame")
plot(New_Vector_r_10)

This is a shapefile now, but for my purpose its fine. One can always convert this back to raster as well.

kjtheron
  • 191
  • 1
  • 11
  • One can edit the ```snap``` parameter in the ```poly2nb``` to select neighbors at appropriate distances. – kjtheron Mar 06 '21 at 05:02