Share

7.2 Landscape metrics within polygons

Some research designs may need landscape metrics for several areas that may be available as a shapefile or some other polygon layer. The code that follows enables researchers to import a shapefile, extract individual polygons from the shapefile, and clip and generate landscape metrics for each polgyon. Again, it is important to be sure the number of habitat patches, patch sizes, and size of the defined poygons will permit estimation of landscape metrics.

 

  1. Exercise 7.1_7.2- Download and extract zip folder into your preferred location
  2. Set working directory to the extracted folder in R under File - Change dir...
  3. First we need to load the packages needed for the exercise


    library(SDMTools)
    library(raster)
    library(rgdal)
    library(maptools)
    library(sp)
    library(adehabitatHR)
    library(rgeos)
    library(plyr)

  4. Now continuing along within the script "Frag_Auto.R", run code below metrics within a single area
    ### load raster file into R
    raster <- raster("county_hab")
    raster
    ### load PA shapefile into R
    HareCounties <- readOGR(dsn=".", layer="Hare_Counties")
    HareCounties
    plot(HareCounties)
    proj4string(HareCounties)
    proj4string(raster)
    image(raster)
    plot(HareCounties, add=T)
    #Let’s project Counties to habitat just to be safe
    new.crs <-CRS("+proj=lcc +lat_1=41.95 +lat_2=40.88333333333333
    +lat_0=40.16666666666666 + lon_0=-77.75 +x_0=600000 +y_0=0
    +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs")
    county <- spTransform(HareCounties, CRS=new.crs)
    HareCounties <- county
    proj4string(HareCounties)
    proj4string(raster)
    #Matching projections successful!
    row.names(HareCounties)<-as.character(HareCounties$COUNTY_NAM)
    names.polygons<-sapply(HareCounties@polygons, function(x) slot(x,"ID"))
    text(coordinates(HareCounties), labels=sapply(slot(HareCounties, "polygons"), function(i) slot(i, "ID")), cex=0.8)

  5. Now we want to export individual counties by County name (i.e., COUNTY_NAM) as individual shapefiles

    indata <- HareCounties
    innames <- unique(HareCounties@data$COUNTY_NAM)
    innames <- innames[1:2]
    outnames <- innames

    # begin loop to create separate county shapefiles
    for (i in 1:length(innames)){
    data <- indata[which(indata$COUNTY_NAM==innames[i]),]
    if(dim(data)[1] != 0){
    writePolyShape(data,fn=paste(outnames[i],sep="/"),factor2char=TRUE)
    write.table(innames, "List.txt", eol=".shp\n", col.names=FALSE, quote=FALSE, row.names=FALSE)
    }}

  6. We then need to read a list of the counties we want to use in our analysis. Multiple text files with shapefile names can be saved if you want to run separately across states or study regions.

    #Read in a list of shapefiles files from above
    Listshps<-read.table("List.txt",sep="\t",header=F)
    Listshps

  7. Now we use the plyr package to create a single function that runs multiple functions in a single statement.

    shape <- function(Listshps) {
    file <- as.character(Listshps[1,])
    shape <-readShapeSpatial(file)
    mask <- mask(raster,shape)
    ### Calculate the Class statistics in each county
    cl.data <- ClassStat(mask)
    }

  8. The line below also uses the plyr package to run the newly created function (shape) over the list of shapefiles called in with the List statement (Listshps) by each ID ([.(id)]) in the Listshps

    results <- ddply(Listshps, .(id), shape)
    write.table(results, "FragCounty.txt")