Area Design
Tom Kincaid
2020-09-17
Preliminaries
This document presents example GRTS survey designs for an area resource. The area resource used in the designs is Omernik level 3 ecoregions within Utah. Four survey designs will be presented: (1) an unstratified, equal probability design; (2) an unstratified, unequal probability design; (3) a stratified, equal probability design; and (4) an unstratified, unequal probability design with an oversample and a panel structure for survey over time. The sampling frame used for the survey designs is contained in either an ESRI shapefile, an sf package object, or an sp package object. The frame contains the coordinates for a set of polygons that define the area resource in addition to attribute data associated with the polygons. The coordinate system for the set of points in the sampling frame is an equal area projection rather than latitude and longitude. An equal area projection is used so that calculation of distance between points is valid.
The initial step is to use the library function to load the spsurvey package. After the package is loaded, a message is printed to the R console indicating that the spsurvey package was loaded successfully.
Load the spsurvey package:
library(spsurvey)
library(sf)
library(sp)
Read the sf object
For creating a survey design using the spsurvey package, the standard form of input regarding the resource is a simple features (sf) object. An sf data set for creating the survey designs in this vignette is included in the data directory of the package. The data function is used to load the data set stored in the data directory into an object named UT_ecoregions. Note that sf objects loaded from the data sets in the data directory are stored in a format that is defined in the sf package. See documentation for the sf package for additional information regarding format of those objects.
Attribute data
The ecoregion attribute will be used to define stratum codes and unequal selection probability (multidensity) categories for the survey designs. Ecoregion is contained in a variable named “Level3_Nam” and includes seven unique values. Frame area is summarized for the ecoregion attribute. Note that ecoregion area measured in hectares is contained in the variable named “Area_ha”. The tapply function is used to calculate total area for each ecoregion. The addmargins function is applied to the output from tapply to calculate total area for all ecoregions, and the round function is used to round value to whole numbers. Finally, the resulting table is displayed.
Display the initial six features in the sf object:
head(UT_ecoregions)
#> Simple feature collection with 6 features and 3 fields
#> geometry type: MULTIPOLYGON
#> dimension: XY
#> bbox: xmin: -1574643 ymin: 1629173 xmax: -1085382 ymax: 2250764
#> projected CRS: NAD83 / Conus Albers
#> Level3 Level3_Nam Area_ha geometry
#> 1 80 Northern Basin and Range 1.42202e+11 MULTIPOLYGON (((-1317868 22...
#> 2 18 Wyoming Basin 1.33312e+11 MULTIPOLYGON (((-1128931 20...
#> 3 13 Central Basin and Range 3.09949e+11 MULTIPOLYGON (((-1487359 21...
#> 4 19 Wasatch and Uinta Mountains 4.47240e+10 MULTIPOLYGON (((-1288590 22...
#> 5 20 Colorado Plateaus 1.26379e+11 MULTIPOLYGON (((-1087925 20...
#> 6 21 Southern Rockies 5.40909e+08 MULTIPOLYGON (((-1133832 18...
Summarize frame area by ecoregion:
temp <- with(UT_ecoregions, tapply(Area_ha, Level3_Nam, sum))
temp <- round(addmargins(temp), 0)
temp
#> Central Basin and Range Colorado Plateaus
#> 309949000000 126379000000
#> Mojave Basin and Range Northern Basin and Range
#> 129599000000 142202000000
#> Southern Rockies Wasatch and Uinta Mountains
#> 946443000 45693759000
#> Wyoming Basin Sum
#> 133312000000 888081202000
The seven ecoregions in Utah are displayed below. Ecoregions are used during creation of several of the survey designs described in this vignette.
Unstratified, equal probability, GRTS survey design
The first survey design is an unstratified, equal probability design. The set.seed function is called so that, if necessary, the designs can be replicated.
The initial step is to create a list named Equaldsgn that contains information for specifying the survey design. Since the survey design is unstratified, the list contains a single item named “None” that also is a list. The “None” list includes two items: panel, which is used to specify the sample size for each panel, and seltype, which is used to input the type of random selection for the design. For this example, panel is assigned a single value named “PanelOne” that is set equal to 115, and seltype is assigned the value “Equal”, which indicates equal probability selection.
The grts function in the spsurvey package is called to select the survey design. The following arguments are included in the call to grts: (1) design: the named list of stratum design specifications, which is assigned the Equaldsgn list; (2) DesignID: name for the design, which is used to create a site ID for each site and is assigned the value “EQUAL”; (3) type.frame: the type of frame, which is assigned the value “area” to indicate an area resource; (4) src.frame: source of the frame, which is assigned the value “sf.object” to indicate an sf object frame; (5) sf.object: the sf object, which is assigned the value UT_ecoregions; and (6) shapefile: option to create a shapefile containing the survey design information, which is assigned FALSE.
During execution of the grts function, messages are printed that indicate the initial number of hierarchical levels used for the GRTS grid, the current number of levels, and the final number of levels. The set of messages is printed for each stratum, and is labeled with the stratum name. For this example, the set of messages is labeled “None”, i.e., the name used in the Equaldsgn list. Upon completion of the call to grts, the initial six sites for the survey design and a design summary are printed. The output object created by the grts function is assigned class “SpatialDesign”. The design summary is created using the summary method for that class. In addition to summary, a plot method is available for the SpatialDesign class. For assistance using the summary and plot methods, see documentation for “SpatialDesign-class” on the R help page for spsurvey.
Call the set.seed function so that the design can be replicated:
Create the design list:
Equaldsgn <- list(None=list(panel=c(PanelOne=50), seltype="Equal"))
Select the sample:
Equalsites <- grts(design=Equaldsgn,
DesignID="EQUAL",
type.frame="area",
src.frame="sf.object",
sf.object=UT_ecoregions,
maxlev = 5,
shapefile=FALSE)
#>
#> Stratum: None
#> Current number of levels: 3
#> Current number of levels: 4
#> Final number of levels: 4
Print the initial six lines of the survey design:
head(Equalsites)
#> coordinates siteID xcoord ycoord mdcaty wgt stratum
#> 1 (-1155597, 2068384) EQUAL-01 -1155597 2068384 Equal 4396596870 None
#> 2 (-1425210, 2065981) EQUAL-02 -1425210 2065981 Equal 4396596870 None
#> 3 (-1198925, 1660384) EQUAL-03 -1198925 1660384 Equal 4396596870 None
#> 4 (-1489064, 1970077) EQUAL-04 -1489064 1970077 Equal 4396596870 None
#> 5 (-1118314, 2031373) EQUAL-05 -1118314 2031373 Equal 4396596870 None
#> 6 (-1380084, 1984276) EQUAL-06 -1380084 1984276 Equal 4396596870 None
#> panel EvalStatus EvalReason Level3 Level3_Nam Area_ha
#> 1 PanelOne NotEval 19 Wasatch and Uinta Mountains 4.47240e+10
#> 2 PanelOne NotEval 13 Central Basin and Range 3.09949e+11
#> 3 PanelOne NotEval 20 Colorado Plateaus 1.26379e+11
#> 4 PanelOne NotEval 13 Central Basin and Range 3.09949e+11
#> 5 PanelOne NotEval 20 Colorado Plateaus 1.26379e+11
#> 6 PanelOne NotEval 13 Central Basin and Range 3.09949e+11
Print the survey design summary:
summary(Equalsites)
#>
#>
#> Design Summary: Number of Sites
#>
#> stratum
#> None Sum
#> 50 50
Unstratified, unequal probability, GRTS survey design
The second survey design is an unstratified, unequal probability design. Ecoregions are used to identify multidensity categories. List Unequaldsgn is assigned design specifications. Since the survey design is unstratified, Unequaldsgn includes a single list named “None” that contains three items: panel, seltype, and caty.n. The value for panel is the same as for the equal probability design, and seltype is assigned “Unequal” to indicate unequal selection probabilities. The third item, caty.n, assigns sample sizes for each of seven multidensity categories, where ecoregion names are used as the categories. Note that the sum of sample sizes provided in caty.n must equal the value in panel.
For this survey design, a shapefile will be used as the sampling frame. The sf package function st_write is used to create the shapefile. The following arguments are included in the call to grts: (1) design: assigned the Unequaldsgn list; (2) DesignID: assigned the value “UNEQUAL”; (3) type.frame: assigned the value “area”; (4) src.frame: assigned the value “shapefile”; (5) in.shape: assigned the value “UT_ecoregions.shp”; (6) mdcaty: name of the column in the attributes data frame that identifies the unequal probability category for each element in the frame, which is assigned the value “level3_nam”; and (7) shapefile: assigned the value FALSE. Upon completion of the call to grts, the initial six sites for the survey design and a design summary are printed.
st_write(UT_ecoregions, "UT_ecoregions.shp", quiet = TRUE, delete_dsn = TRUE)
#> Warning in CPL_write_ogr(obj, dsn, layer, driver,
#> as.character(dataset_options), : GDAL Error 1: UT_ecoregions.shp does not appear
#> to be a file or directory.
#> Warning in CPL_write_ogr(obj, dsn, layer, driver,
#> as.character(dataset_options), : GDAL Message 1: Value 142202000000 of field
#> Area_ha of feature 0 not successfully written. Possibly due to too larger number
#> with respect to field width
#> Warning in CPL_write_ogr(obj, dsn, layer, driver,
#> as.character(dataset_options), : GDAL Message 1: Value 133312000000 of field
#> Area_ha of feature 1 not successfully written. Possibly due to too larger number
#> with respect to field width
#> Warning in CPL_write_ogr(obj, dsn, layer, driver,
#> as.character(dataset_options), : GDAL Message 1: Value 309949000000 of field
#> Area_ha of feature 2 not successfully written. Possibly due to too larger number
#> with respect to field width
#> Warning in CPL_write_ogr(obj, dsn, layer, driver,
#> as.character(dataset_options), : GDAL Message 1: Value 44724000000 of field
#> Area_ha of feature 3 not successfully written. Possibly due to too larger number
#> with respect to field width
#> Warning in CPL_write_ogr(obj, dsn, layer, driver,
#> as.character(dataset_options), : GDAL Message 1: Value 126379000000 of field
#> Area_ha of feature 4 not successfully written. Possibly due to too larger number
#> with respect to field width
#> Warning in CPL_write_ogr(obj, dsn, layer, driver,
#> as.character(dataset_options), : GDAL Message 1: Value 540909000 of field
#> Area_ha of feature 5 not successfully written. Possibly due to too larger number
#> with respect to field width
#> Warning in CPL_write_ogr(obj, dsn, layer, driver,
#> as.character(dataset_options), : GDAL Message 1: Value 129599000000 of field
#> Area_ha of feature 6 not successfully written. Possibly due to too larger number
#> with respect to field width
#> Warning in CPL_write_ogr(obj, dsn, layer, driver,
#> as.character(dataset_options), : GDAL Message 1: Value 343657000 of field
#> Area_ha of feature 7 not successfully written. Possibly due to too larger number
#> with respect to field width
#> Warning in CPL_write_ogr(obj, dsn, layer, driver,
#> as.character(dataset_options), : GDAL Message 1: Value 626102000 of field
#> Area_ha of feature 8 not successfully written. Possibly due to too larger number
#> with respect to field width
#> Warning in CPL_write_ogr(obj, dsn, layer, driver,
#> as.character(dataset_options), : GDAL Message 1: Value 405534000 of field
#> Area_ha of feature 9 not successfully written. Possibly due to too larger number
#> with respect to field width
Create the design list:
Unequaldsgn <- list(None=list(panel=c(PanelOne=50),
seltype="Unequal",
caty.n=c("Central Basin and Range"=10,
"Colorado Plateaus"=10,
"Mojave Basin and Range"=5,
"Northern Basin and Range"=5,
"Southern Rockies"=5,
"Wasatch and Uinta Mountains"=10,
"Wyoming Basin"=5)))
Select the sample:
Unequalsites <- grts(design=Unequaldsgn,
DesignID="UNEQUAL",
type.frame="area",
src.frame="shapefile",
in.shape="UT_ecoregions.shp",
mdcaty="Level3_Nam",
maxlev = 3,
shapefile=FALSE)
#>
#> Stratum: None
#> Current number of levels: 3
#> Final number of levels: 3
#> Warning in grtsarea(sframe, sum(n.desired), SiteBegin, shift.grid, startlev, :
#> Of the 33 grid cells from which sample points were selected,
#> 12 (36.4%) of the cells contained more than one sample point.
#> Warning in showSRID(uprojargs, format = "PROJ", multiline = "NO"): Discarded
#> datum Unknown based on GRS80 ellipsoid in CRS definition
Print the initial six lines of the survey design:
head(Unequalsites)
#> coordinates siteID xcoord ycoord mdcaty
#> 1 (-1558089, 1703408) UNEQUAL-01 -1558089 1703408 Mojave Basin and Range
#> 2 (-1126103, 2074736) UNEQUAL-02 -1126103 2074736 Wyoming Basin
#> 3 (-1263654, 2192250) UNEQUAL-03 -1263654 2192250 Wyoming Basin
#> 4 (-1162018, 1723076) UNEQUAL-04 -1162018 1723076 Southern Rockies
#> 5 (-1559142, 1724424) UNEQUAL-05 -1559142 1724424 Mojave Basin and Range
#> 6 (-1354725, 2014889) UNEQUAL-06 -1354725 2014889 Central Basin and Range
#> wgt stratum panel EvalStatus EvalReason Level3 Area_ha
#> 1 386372151 None PanelOne NotEval 14 1.29599e+11
#> 2 582117660 None PanelOne NotEval 18 1.33312e+11
#> 3 582117660 None PanelOne NotEval 18 1.33312e+11
#> 4 189288683 None PanelOne NotEval 21 4.05534e+08
#> 5 386372151 None PanelOne NotEval 14 1.29599e+11
#> 6 8206454636 None PanelOne NotEval 13 3.09949e+11
Print the survey design summary:
summary(Unequalsites)
#> Design Summary: Number of Sites Classified by mdcaty (Multidensity Category)
#>
#> mdcaty
#> Central Basin and Range Colorado Plateaus
#> 10 10
#> Mojave Basin and Range Northern Basin and Range
#> 5 4
#> Southern Rockies Wasatch and Uinta Mountains
#> 6 9
#> Wyoming Basin Sum
#> 6 50
Stratified, equal probability, GRTS survey design
The third survey design is a stratified, equal probability design. Ecoregions are used to identify strata. List Stratdsgn is assigned design specifications. The ecoregion attribute is used to identify strata. Stratdsgn includes seven lists, one for each stratum. The names for the lists match the levels of the stratum variable, i.e., the unique values of the ecoregion attribute. Each list in Stratdsgn contains two items: panel and seltype. The value for panel is the same as for the equal probability design, and seltype is assigned “Equal”.
For this survey design, an sp package object will be used as the sampling frame. The sf package function as_Spatial is used to create an sp object named UT_ecoregions_sp. The following arguments are included in the call to grts: (1) design: assigned the Stratdsgn list; (2) DesignID: assigned the value “STRATIFIED”; (3) type.frame: assigned the value “area”; (4) src.frame: assigned the value “sp.object” to indicate that the sampling frame is provided by an sp object; (5) sp.object: the sp object, which is assigned the UT_ecoregions_sp object; (6) stratum: name of the column in the attributes data frame that identifies the stratum code for each element in the frame, which is assigned the value “level3_nam”; and (7) shapefile: assigned the value FALSE. Upon completion of the call to grts, the initial six sites for the survey design and a design summary are printed.
Create the sp object:
UT_ecoregions_sp <- sf::as_Spatial(UT_ecoregions)
#> Warning: st_crs<- : replacing crs does not reproject data; use st_transform for
#> that
proj4string(UT_ecoregions_sp) <- sp::CRS(st_crs(UT_ecoregions)$proj4string)@projargs
#> Warning in proj4string(obj): CRS object has comment, which is lost in output
Create the design list:
Stratdsgn <- list("Central Basin and Range"=list(panel=c(PanelOne=10),
seltype="Equal"),
"Colorado Plateaus"=list(panel=c(PanelOne=10),
seltype="Equal"),
"Mojave Basin and Range"=list(panel=c(PanelOne=5),
seltype="Equal"),
"Northern Basin and Range"=list(panel=c(PanelOne=5),
seltype="Equal"),
"Southern Rockies"=list(panel=c(PanelOne=5),
seltype="Equal"),
"Wasatch and Uinta Mountains"=list(panel=c(PanelOne=10),
seltype="Equal"),
"Wyoming Basin"=list(panel=c(PanelOne=5),
seltype="Equal"))
Select the sample:
Stratsites <- grts(design=Stratdsgn,
DesignID="STRATIFIED",
type.frame="area",
src.frame="sp.object",
sp.object=UT_ecoregions_sp,
stratum="Level3_Nam",
maxlev = 3,
shapefile=FALSE)
#>
#> Stratum: Central Basin and Range
#> Current number of levels: 2
#> Current number of levels: 3
#> Final number of levels: 3
#>
#> Stratum: Colorado Plateaus
#> Current number of levels: 2
#> Current number of levels: 3
#> Final number of levels: 3
#>
#> Stratum: Mojave Basin and Range
#> Current number of levels: 2
#> Final number of levels: 2
#>
#> Stratum: Northern Basin and Range
#> Current number of levels: 2
#> Current number of levels: 3
#> Final number of levels: 3
#>
#> Stratum: Southern Rockies
#> Current number of levels: 2
#> Current number of levels: 3
#> Final number of levels: 3
#>
#> Stratum: Wasatch and Uinta Mountains
#> Current number of levels: 2
#> Current number of levels: 3
#> Final number of levels: 3
#>
#> Stratum: Wyoming Basin
#> Current number of levels: 2
#> Current number of levels: 3
#> Final number of levels: 3
Print the initial six lines of the survey design:
head(Stratsites)
#> coordinates siteID xcoord ycoord mdcaty wgt
#> 1 (-1531721, 1824033) STRATIFIED-01 -1531721 1824033 Equal 8206454636
#> 2 (-1398588, 2056297) STRATIFIED-02 -1398588 2056297 Equal 8206454636
#> 3 (-1369118, 2007008) STRATIFIED-03 -1369118 2007008 Equal 8206454636
#> 4 (-1338733, 1981828) STRATIFIED-04 -1338733 1981828 Equal 8206454636
#> 5 (-1443186, 2136787) STRATIFIED-05 -1443186 2136787 Equal 8206454636
#> 6 (-1353950, 2214941) STRATIFIED-06 -1353950 2214941 Equal 8206454636
#> stratum panel EvalStatus EvalReason Level3 Area_ha
#> 1 Central Basin and Range PanelOne NotEval 13 3.09949e+11
#> 2 Central Basin and Range PanelOne NotEval 13 3.09949e+11
#> 3 Central Basin and Range PanelOne NotEval 13 3.09949e+11
#> 4 Central Basin and Range PanelOne NotEval 13 3.09949e+11
#> 5 Central Basin and Range PanelOne NotEval 13 3.09949e+11
#> 6 Central Basin and Range PanelOne NotEval 13 3.09949e+11
Print the survey design summary:
summary(Stratsites)
#>
#>
#> Design Summary: Number of Sites
#>
#> stratum
#> Central Basin and Range Colorado Plateaus
#> 10 10
#> Mojave Basin and Range Northern Basin and Range
#> 5 5
#> Southern Rockies Wasatch and Uinta Mountains
#> 5 10
#> Wyoming Basin Sum
#> 5 50
Unstratified, unequal probability, GRTS survey design with an oversample and a panel structure for survey over time
The fourth survey design is an unstratified, unequal probability design with an oversample and a panel structure for survey over time. List Paneldsgn is assigned design specifications. Since the survey design is unstratified, Paneldsgn includes a single list named “None” that contains four items: panel, seltype, caty.n, and over. A vector identifying sample sizes for five panels is assigned to panel. The value “Unequal” is assigned to seltype, which indicates unequal selection probabilities. The third item, caty.n, assigns sample sizes for each of seven multidensity categories, where ecoregion names are used as the categories. Note that the sum of sample sizes provided in caty.n must equal the sum of sample sizes in panel. The value 100 is assigned to over, which specifies an oversample of 100 sites. An oversample is replacement sites for the survey design. The grts function attempts to distribute the oversample proportionately among sample sizes for the multidensity categories. If the oversample proportion for one or more categories is not a whole number, a warning message is printed and the proportion is rounded to the next higher integer. For this example, the oversample is not proportionate to the category sample sizes, and the warning message is printed by calling the warnings function.
For this survey design, an sf object will be used as the sampling frame. The following arguments are included in the call to grts: (1) design: assigned the Paneldsgn list; (2) DesignID: assigned the value “UNEQUAL”; (3) type.frame: assigned the value “area”; (4) src.frame: assigned the value “sf.object”; (5) sf.object: the sf object, which is assigned the value UT_ecoregions; (6) mdcaty: assigned the value “level3_nam”; and (7) shapefile: assigned the value FALSE. Upon completion of the call to grts, the initial six sites for the survey design and a design summary are printed.
Create the design list:
Paneldsgn <- list(None=list(panel=c(Year1=10, Year2=10, Year3=10,
Year4=10, Year5=10),
seltype="Unequal",
caty.n=c("Central Basin and Range"=10,
"Colorado Plateaus"=10,
"Mojave Basin and Range"=5,
"Northern Basin and Range"=5,
"Southern Rockies"=5,
"Wasatch and Uinta Mountains"=10,
"Wyoming Basin"=5),
over=5))
Select the sample:
Panelsites <- grts(design=Paneldsgn,
DesignID="UNEQUAL",
type.frame="area",
src.frame="sf.object",
sf.object=UT_ecoregions,
maxlev = 5,
mdcaty="Level3_Nam",
shapefile=FALSE)
#>
#> Stratum: None
#> Warning in grts(design = Paneldsgn, DesignID = "UNEQUAL", type.frame = "area", :
#> Oversample size is not proportional to category sample sizes for stratum
#> "None".
#> Current number of levels: 3
#> Current number of levels: 5
#> Final number of levels: 5
#> Warning in grtsarea(sframe, sum(n.desired), SiteBegin, shift.grid, startlev, :
#> Of the 55 grid cells from which sample points were selected,
#> 2 (3.6%) of the cells contained more than one sample point.
Print the warning message:
Print the initial six lines of the survey design:
head(Panelsites)
#> coordinates siteID xcoord ycoord mdcaty
#> 1 (-1360080, 2099839) UNEQUAL-01 -1360080 2099839 Central Basin and Range
#> 2 (-1120360, 2076983) UNEQUAL-02 -1120360 2076983 Wyoming Basin
#> 3 (-1144710, 1796206) UNEQUAL-03 -1144710 1796206 Southern Rockies
#> 4 (-1141340, 1805220) UNEQUAL-04 -1141340 1805220 Southern Rockies
#> 5 (-1249382, 2133024) UNEQUAL-05 -1249382 2133024 Wyoming Basin
#> 6 (-1383508, 2067069) UNEQUAL-06 -1383508 2067069 Central Basin and Range
#> wgt stratum panel EvalStatus EvalReason Level3 Area_ha
#> 1 8206454636 None Year1 NotEval 13 3.09949e+11
#> 2 582117660 None Year1 NotEval 18 1.33312e+11
#> 3 189288683 None Year1 NotEval 21 5.40909e+08
#> 4 189288683 None Year1 NotEval 21 5.40909e+08
#> 5 582117660 None Year1 NotEval 18 1.33312e+11
#> 6 8206454636 None Year1 NotEval 13 3.09949e+11
Print the survey design summary:
summary(Panelsites)
#>
#>
#> Design Summary: Number of Sites Classified by mdcaty (Multidensity Category)
#> and panel
#>
#> panel
#> mdcaty OverSamp Year1 Year2 Year3 Year4 Year5 Sum
#> Central Basin and Range 1 3 4 0 2 2 12
#> Colorado Plateaus 3 1 1 2 3 3 13
#> Mojave Basin and Range 1 1 0 2 0 1 5
#> Northern Basin and Range 0 1 1 2 2 1 7
#> Southern Rockies 1 2 1 0 1 1 6
#> Wasatch and Uinta Mountains 0 0 3 3 2 2 10
#> Wyoming Basin 1 2 0 1 0 0 4
#> Sum 7 10 10 10 10 10 57
