spatstat ppp objects from FIA?

[jgrn307]

I'm intending to use ppp objects for spatial point pattern analysis, and one of the requirements is to be able to associate the observation window with the individual tree lists. Is there a current way to convert a set of FIA plots into ppp objects? This requires the tree x,y plus marks as well as the plot boundary as a polygon (I believe they are circular plots, right?)

[hunter-stanke]

Sorry for the (very) delayed reply here! I accidentally shut off my email notifications, and I'm just now seeing issues that were opened in August.

I haven't worked with ppp objects before, however I do have some code that will produce spatial polygons of FIA plot boundaries, and I will follow up with some code that will produce a point pattern (sf) of tree locations within plots.

For the plot boundaries, how does this look?

library(rFIA)
library(dplyr)
library(sf)

## Using the canned RI subset
data(fiaRI)

## Get spatial information for all of most recent plot visits
plts <- fiaRI %>% 
  clipFIA() %>%
  area(byPlot = TRUE, returnSpatial = TRUE) %>%
  ## Albers equal area (meters)
  st_transform(crs = 5070)

## Get coordinates of center subplot as a dataframe
c.coords <- st_coordinates(plts) %>%
  as.data.frame() %>%
  ## Label the center subplot appropriately
  mutate(plot = 1:n(), 
         subplot = 1)

## Imagine the non-center subplot coordinates as the vertices of an 
## equilateral triangle, where sides are of length sqrt(120^2 + 120^2) dt
## We can then compute the x and y coordinates of non-center subplots directly
side.length = sqrt(120^2 + 120^2) * 0.3048 # meters

## Compute coordinates of subplot 2
t.coords <- c.coords %>%
  ## Shift X 120ft up
  mutate(Y = Y + (120 * 0.3048)) %>%
  mutate(subplot = 2)

## Coordinates of bottom subplots
bl.coords <- c.coords %>%
  ## Some trig
  mutate(Y = Y - (side.length / 2),
         X = X - (side.length / 2)) %>% 
  mutate(subplot = 3)
br.coords <- bl.coords %>%
  ## Some trig
  mutate(X = X + side.length) %>% 
  mutate(subplot = 4)

## Combine all coordinates
coords <- bind_rows(c.coords, t.coords, bl.coords, br.coords) 

## Use spatial buffers to make polygons delineating micro-, sub-, and macro-plots
# subplots
subp <- coords %>%
  st_as_sf(coords = c('X', 'Y')) %>%
  st_buffer(dist = 24 * 0.3048)

# microplots
micr <- coords %>%
  mutate(X = X + (12 * 0.3048)) %>%
  st_as_sf(coords = c('X', 'Y')) %>%
  st_buffer(dist = 6.8 * 0.3048)

# macroplots
macr <- coords %>%
  st_as_sf(coords = c('X', 'Y')) %>%
  st_buffer(dist = 58.9 * 0.3048)

## Set projection appropriate prior to writing
st_crs(subp) <- 5070
st_crs(micr) <- 5070
st_crs(macr) <- 5070

## Make sure things look right
library(ggplot2)
ggplot() +
  geom_sf(data = filter(macr, plot == 1)) +
  geom_sf(data = filter(subp, plot == 1)) +
  geom_sf(data = filter(micr, plot == 1))

[hunter-stanke]

Follow up from my previous, how does this look for the tree-level point patterns?

library(rFIA)
library(dplyr)
library(sf)

## Using the canned RI subset
data(fiaRI)

## Get spatial information for all of most recent plot visits
plts <- fiaRI %>% 
  clipFIA() %>%
  area(byPlot = TRUE, returnSpatial = TRUE) %>%
  ## Albers equal area (meters)
  st_transform(crs = 5070)

## Get coordinates of center subplot as a dataframe
c.coords <- st_coordinates(plts) %>%
  as.data.frame() %>%
  ## Label the center SUBP appropriately
  mutate(pltID = plts$pltID, 
         SUBP = 1)

## Imagine the non-center SUBP coordinates as the vertices of an 
## equilateral triangle, where sides are of length sqrt(120^2 + 120^2) dt
## We can then compute the x and y coordinates of non-center SUBPs directly
side.length = sqrt(120^2 + 120^2) * 0.3048 # meters

## Compute coordinates of SUBP 2
t.coords <- c.coords %>%
  ## Shift X 120ft up
  mutate(Y = Y + (120 * 0.3048)) %>%
  mutate(SUBP = 2)

## Coordinates of bottom SUBPs
bl.coords <- c.coords %>%
  ## Some trig
  mutate(Y = Y - (side.length / 2),
         X = X - (side.length / 2)) %>% 
  mutate(SUBP = 3)
br.coords <- bl.coords %>%
  ## Some trig
  mutate(X = X + side.length) %>% 
  mutate(SUBP = 4)

## Combine all coordinates
coords <- bind_rows(c.coords, t.coords, bl.coords, br.coords) %>%
  rename(X.PLOT = X, Y.PLOT = Y)

## Use spatial buffers to make polygons subplots
subp <- coords %>%
  st_as_sf(coords = c('X.PLOT', 'Y.PLOT')) %>%
  st_buffer(dist = 24 * 0.3048)

## Now get the offsets in tree location from plot center
tree <- fiaRI %>% 
  clipFIA() %>%
  tpa(grpBy = c(SUBP, TREE, DIST, AZIMUTH), byPlot = TRUE) %>%
  ## Subplot only
  filter(TPA == 6.018046) %>%
  ## Convert degrees to radians
  ## Convert feet to meters
  mutate(rad = AZIMUTH * (pi/180),
         DIST = DIST * 0.3048) %>%
  mutate(x = case_when(AZIMUTH %in% c(0, 180) ~ 0,
                       AZIMUTH %in% c(90, 270) ~ DIST,
                       AZIMUTH < 90 ~ sin(rad) * DIST,
                       AZIMUTH < 180 ~ sin(pi - rad) * DIST,
                       AZIMUTH < 270 ~ -sin(rad - pi) * DIST,
                       AZIMUTH < 360 ~ -sin(2*pi - rad) * DIST),
         y = case_when(AZIMUTH %in% c(0, 180) ~ 0,
                       AZIMUTH %in% c(90, 270) ~ DIST,
                       AZIMUTH < 90 ~ cos(rad) * DIST,
                       AZIMUTH < 180 ~ -cos(pi - rad) * DIST,
                       AZIMUTH < 270 ~ -cos(rad - pi) * DIST,
                       AZIMUTH < 360 ~ cos(2*pi - rad) * DIST)) %>%
  select(pltID, 
         SUBP,
         X.TREE = x, 
         Y.TREE = y)

## Make tree coordinates absolute and make spatial
tree.coords <- coords %>%
  left_join(tree, by = c('pltID', 'SUBP')) %>%
  mutate(X = X.PLOT + X.TREE,
         Y = Y.PLOT + Y.TREE) %>%
  ## Drop non-treed forested plots
  filter(!is.na(X) | !is.na(Y)) %>%
  st_as_sf(coords = c('X', 'Y'))

## Check it out for a single plot
library(ggplot2)
ggplot() +
  geom_sf(data = filter(subp, pltID == '1_44_1_228')) +
  geom_sf(data = filter(tree.coords, pltID == '1_44_1_228')) 

[jgrn307]

Hunter, thanks a ton for looking into this -- stay tuned, I think we developed (independently) a workflow working and after I confirm, I can tweak into a functional form and send to you as a potential add-on to rFIA if you'd like!

[hunter-stanke]

That would be excellent - we're always open to extensions!