I have a shapefile with several MULTILINESTRING geometries, representing the borders of some species' distributions, and I want to count the number the number of species occurring in each cell.

I used the rasterize function from the terra package:

count.sp <- terra::rasterize(sp.lines, rast.w, fun = "count")

As an example, the lines in the figure below represent the borders of the distribution of a given species "X":

Example with a species whose distribution is subdivided into several elements

In the example, species' distribution was subdivided into several disjoint parts crossing the same raster pixel. With the code above, I got as a result the number of parts in each pixel (8 and 5, respectively), whereas I would like to get all the parts belonging to the same species and crossing the same pixel be counted as a single occurrence (i.e., I would like to get 1 in both pixels in the example below).

I also tried the function "sum" (with no field, and I obtained 208 and 130, I don't understand what these values represent)

You can find some example data here (1 shapefile with just 2 species and the raster)

What I am missing?

  • What answer do you want given the above data? I don't understand how you get "circles" from lines. In your shapefile, each row is a MULTILINESTRING object. Does one row = one species, and because its MULTI can be made of lots of little disjoint parts. Do you want to add 1 for every grid cell than any part of a MULTI-part object occurs in?
    – Spacedman
    May 9 at 14:25
  • 1
    I would take a look at terra::extract rather than trying to rasterize. With the cells=TRUE argument, unique cell ID's are also returned which, allows for an aggregated count. May 9 at 14:27
  • @Spacedman, "circles" are the borders of species' distributions, obtained from of original polygons after converting them to multilinestring elements. I am using multilines geometries instead of original polygons because I want to count species occurrences even if the distribution does not touch the center of the pixel.
    – M.S.
    May 9 at 14:37
  • @Spacedman, I would like all disjoint parts belonging to the same species within a same pixel be counted as a single occurrence. In other words, for the distribution of the species in the figure above, I would like to get 1 in both pixels, insted of the number of disjoint parts in each pixel (that are 8 on the left and 5 on the right).
    – M.S.
    May 9 at 14:38
  • @JeffreyEvans, thank you for the suggestion, I achieved it following your hint, even if it would be simpler if it was possible directly with rasterize
    – M.S.
    May 9 at 16:09

1 Answer 1


Following @JeffreyEvans suggestion, I achieved it using an "extract" approach:

First, I assigned an unique ID to each raster pixel:

values(rast.w) <- 1:ncell(rast.w)

Then, I used the extract function from the package terra to get the ID of the cells crossed by each species:

sp <- terra::extract(rast.w, sp.lines)

Then, I converted the resulting list into a dataframe with the ID of each cell (first column) and the number of species occurring in each of them (second column, aggregated counting)

count <- as.data.frame(table(unlist(sp)))

Finally, I created a new raster from the first one, replacing ID values with the counting result:

newraster <- rast.w
newraster <- subs(rast.w, count, by = "test")

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.