[R-sig-eco] How to extract values from contiguous raster cells that ate not touched by SpatialLines?

[Andre Rovai]

Hi all,

I've been trying to extract values from a single attribute raster (area, in
m2) that overlaps with lines (that is, a .shp SpatialLines).

The problem is that, along these lines, my raster sometimes goes from one
to several contiguous cells in all directions.  Using the
extract function only values from cells that are touched by the lines are
extracted.  Thus, when I add up the extracted values from all lines a
significant amount of area (m2) is lost due to cells that were not touched
by the line and therefore values were not extracted.

I tried to work it around by:

Step 1 - first aggregating my raster to a lower resolution (i.e. increasing
the fact argument) and then
Step 2 - rasterizing the lines using this aggregated raster (created in
step 1) as a mold to make sure the rasterized lines would get thick enough
to cover the horizontal spread of cells in my original resolution raster.
Step 3 - Then I resample the rasterized lines (created in step 2) back to
the original resolution I started with.
Step 4 - Finally, extracted the values from the resampled rasterized lines
(created in step 3).

However, it didn't quite work as now the total area (m2) varies according
to the fact="" value I use when first aggregating the raster (in step 1).

I really appreciate if anyone has already dealt with a similar problem and
can help me out here.  Here are the codes I've been running to try to get
it to work:


# input raster file

g.025 <- raster("ras.asc")

g.1 <- aggregate(g.025, fact=2, fun=sum)



# input SpatialLines

Spline1 <- readOGR("/Users/xxxxx.shp")

Spline2 <- readOGR("/Users/xxxxx.shp")

Spline3 <- readOGR("/Users/xxxxx.shp")



# rasterizing using low resolution raster (aggregated)

c1 <- rasterize(Spline1, g.1, field=Spline1$type, fun=sum)

c2 <- rasterize(Spline2, g.1, field=Spline2$type, fun=sum)

c3 <- rasterize(Spline3, g.1, field=Spline3$type, fun=sum)



# resampling back to higher resolution

c1 <- resample(c1, g.025)

c2 <- resample(c2, g.025)

c3 <- resample(c3, g.025)



# preparing to extract area (m2) values from raster “g.025”

c1tab <- as.data.frame(c1, xy=T)

c2tab <- as.data.frame(c2, xy=T)

c3tab <- as.data.frame(c3, xy=T)

c1tab <- c1tab[which(is.na(c1tab$layer)!=T),]

c2tab <- c2tab[which(is.na(c2tab$layer)!=T),]

c3tab <- c3tab[which(is.na(c3tab$layer)!=T),]



# extracting area (m2) values from raster “g.025”

c1tab[,4] <- extract(g.025, c1tab[,1:2])

c2tab[,4] <- extract(g.025, c2tab[,1:2])

c3tab[,4] <- extract(g.025, c3tab[,1:2])

names(c1tab)[4] <- "area_m2"

names(c2tab)[4] <- "area_m2"

names(c3tab)[4] <- "area_m2"



# sum total area (m2)

c1_area <- sum(c1tab$area_m2)

c2_area <- sum(c2tab$area_m2)

c3_area <- sum(c3tab$area_m2)

tot_area <- sum(c1_area, c2_area, c3_area)


Thanks!

Andre Rovai
Department of Oceanography and Coastal Sciences
Louisiana State University

	[[alternative HTML version deleted]]