Difficulty using huc watersheds as polygons
MalakoffE opened this issue · 2 comments
Hello!
I am running into trouble using VT HUC_12_watershed polygons in a tpa function. The code runs, but outputs a table with 26 variables and no observations. The polygons are stored in a SpatialPolygonsDataFrame. Maybe they are too geographically small? My code (which is super straightforward) is below.
set wd
setwd("/Users/Liza/OneDrive - Dartmouth College/BigDataHydrology/Final_Project")
Clear environment
rm(list = ls())
Libraries
library(rFIA)
library(RSQLite)
library(data.table)
library(sf)
library(rgdal)
library(here)
READ FIA DATA
Downloaded from FIA datahub SQL file using fix from Github
And keep only most recent data from each plot
fiaVT = readFIA("/Users/Liza/OneDrive - Dartmouth College/BigDataHydrology/Final_Project/VT_FIA_data", common = TRUE)
fiaVT = fiaVT_recent <- clipFIA(fiaVT, mostRecent = TRUE)
GET WATERSHED POLY
Downloaded from https://geodata.vermont.gov/datasets/vt-subwatershed-boundaries-huc12/explore
huc12_vt <- readOGR(dsn = '/Users/Liza/Library/CloudStorage/OneDrive-DartmouthCollege/BigDataHydrology/Final_Project/Watershed_shapefiles/VT_Subwatershed_Boundaries_-HUC12', layer = 'VT_Subwatershed_Boundaries-_HUC12')
Calculate tpa for each watershed
tpaHUC <- tpa(fiaVT, polys = huc12_vt, returnSpatial = TRUE)
You’re likely right that the HUC12 watersheds are too small, particularly for Vermont with such a small area already. My suggestion is create a spatial points object with the FIA data (tpa( byPlot=T, returnSpatial=T) and overlay those on the HUC12 polygons to see how many FIA plots are in each watershed. Likely not enough to make an aggregate estimate.
More generally, a rule of thumb with FIA data is the smallest unit of aggregation recommended is the county. And in some cases even that is too small depending on tolerance for sample error.