trend_analysis {spsurvey} | R Documentation |
Trend analysis
Description
This function organizes input and output for estimation of trend across time
for a series of samples (for categorical and continuous variables). Trend is estimated using the
analytical procedure identified by the model arguments. For categorical
variables, the choices for the model_cat
argument are: (1) simple linear
regression, (2) weighted linear regression, and (3) generalized linear
mixed-effects model. For continuous variables, the choices for the
model_cont
argument are: (1) simple linear regression, (2) weighted
linear regression, and (3) linear mixed-effects model. The analysis data,
dframe
, can be either a data frame or a simple features (sf
) object. If an
sf
object is used, coordinates are extracted from the geometry column in the
object, arguments xcoord
and ycoord
are assigned values
"xcoord"
and "ycoord"
, respectively, and the geometry column is
dropped from the object.
Usage
trend_analysis(
dframe,
vars_cat = NULL,
vars_cont = NULL,
subpops = NULL,
model_cat = "SLR",
cat_rhs = NULL,
model_cont = "LMM",
cont_rhs = NULL,
siteID = "siteID",
yearID = "year",
weight = "weight",
xcoord = NULL,
ycoord = NULL,
stratumID = NULL,
clusterID = NULL,
weight1 = NULL,
xcoord1 = NULL,
ycoord1 = NULL,
sizeweight = FALSE,
sweight = NULL,
sweight1 = NULL,
fpc = NULL,
popsize = NULL,
invprboot = TRUE,
nboot = 1000,
vartype = "Local",
jointprob = "overton",
conf = 95,
All_Sites = FALSE
)
Arguments
dframe |
Data to be analyzed (analysis data). A data frame or
|
vars_cat |
Vector composed of character values that identify the names
of categorical response variables in |
vars_cont |
Vector composed of character values that identify the
names of continuous response variables in |
subpops |
Vector composed of character values that identify the
names of subpopulation (domain) variables in |
model_cat |
Character value identifying the analytical procedure used
for trend estimation for categorical variables. The choices are:
|
cat_rhs |
Character value specifying the right hand side of the formula
for a generalized linear mixed-effects model. If a value is not provided,
the argument is assigned a value that specifies the Piepho and Ogutu (2002)
model. The default value is |
model_cont |
Character value identifying the analytical procedure used
for trend estimation for continuous variables. The choices are:
|
cont_rhs |
Character value specifying the right hand side of the
formula for a linear mixed-effects model. If a value is not provided, the
argument is assigned a value that specifies the Piepho and Ogutu (2002)
model. The default value is |
siteID |
Character value providing name of the site ID variable in
|
yearID |
Character value providing name of the time period variable in
|
weight |
Character value providing name of the design weight
variable in |
xcoord |
Character value providing name of the x-coordinate variable in
|
ycoord |
Character value providing name of the y-coordinate variable in
|
stratumID |
Character value providing name of the stratum ID variable in
|
clusterID |
Character value providing name of the cluster (stage one) ID
variable in |
weight1 |
Character value providing name of the stage one weight
variable in |
xcoord1 |
Character value providing name of the stage one x-coordinate
variable in |
ycoord1 |
Character value providing name of the stage one y-coordinate
variable in |
sizeweight |
Logical value that indicates whether size weights should be
used during estimation, where |
sweight |
Character value providing name of the size weight variable in
|
sweight1 |
Character value providing name of the stage one size weight
variable in |
fpc |
Object that specifies values required for calculation of the finite population correction factor used during variance estimation. The object must match the survey design in terms of stratification and whether the design is single-stage or two-stage. For an unstratified design, the object is a vector. The vector is composed of a single numeric value for a single-stage design. For a two-stage unstratified design, the object is a named vector containing one more than the number of clusters in the sample, where the first item in the vector specifies the number of clusters in the population and each subsequent item specifies the number of stage two units for the cluster. The name for the first item in the vector is arbitrary. Subsequent names in the vector identify clusters and must match the cluster IDs. For a stratified design, the object is a named list of vectors, where names must match the strata IDs. For each stratum, the format of the vector is identical to the format described for unstratified single-stage and two-stage designs. Note that the finite population correction factor is not used with the local mean variance estimator. Example fpc for a single-stage unstratified survey design:
Example fpc for a single-stage stratified survey design:
Example fpc for a two-stage unstratified survey design:
Example fpc for a two-stage stratified survey design:
|
popsize |
Object that provides values for the population argument of the
Example popsize for calibration:
Example popsize for post-stratification using a data frame:
Example popsize for post-stratification using a table:
Example popsize for post-stratification using an xtabs object:
|
invprboot |
Logical value that indicates whether the inverse probability
bootstrap procedure is used to calculate trend parameter estimates. This
bootstrap procedure is only available for the "LMM" option for continuous
variables. Inverse probability references the design weights, which
are the inverse of the sample inclusion probabilities. The default value
is |
nboot |
Numeric value for the number of bootstrap iterations. The
default is |
vartype |
Character value providing choice of the variance estimator,
where |
jointprob |
Character value providing choice of joint inclusion
probability approximation for use with Horvitz-Thompson and Yates-Grundy
variance estimators, where |
conf |
Numeric value for the Gaussian-based confidence level. The default is
|
All_Sites |
A logical variable used when |
Value
The analysis results. A list composed of two data frames containing trend estimates for all combinations of population Types, subpopulations within Types, and response variables. For categorical variables, trend estimates are calculated for each category of the variable. The two data frames in the output list are:
catsum
data frame containing trend estimates for categorical variables
contsum
data frame containing trend estimates for continuous variables
For the SLR and WLR model options, the data frame contains the following variables:
- Type
subpopulation (domain) name
- Subpopulation
subpopulation name within a domain
- Indicator
response variable
- Trend_Estimate
trend estimate
- Trend_Std_Error
trend standard error
- Trend_LCBxxPct
trend xx% (default 95%) lower confidence bound
- Trend_UCBxxPct
trend xx% (default 95%) upper confidence bound
- Trend_p_Value
trend p-value
- Intercept_Estimate
intercept estimate
- Intercept_Std_Error
intercept standard error
- Intercept_LCBxxPct
intercept xx% (default 95%) lower confidence bound
- Intercept_UCBxxPct
intercept xx% (default 95%) upper confidence bound
- Intercept_p_Value
intercept p-value
- R_Squared
R-squared value
- Adj_R_Squared
adjusted R-squared value
For the GLMM and LMM model options, contents of the data frames will vary
depending on the model specified by arguments cat_rhs
and
cont_rhs
. For the default PO model, the data frame contains the
following variables:
- Type
subpopulation (domain) name
- Subpopulation
subpopulation name within a domain
- Indicator
response variable
- Trend_Estimate
trend estimate
- Trend_Std_Error
trend standard error
- Trend_LCBxxPct
trend xx% (default 95%) lower confidence bound
- Trend_UCBxxPct
trend xx% (default 95%) upper confidence bound
- Trend_p_Value
trend p-value
- Intercept_Estimate
intercept estimate
- Intercept_Std_Error
intercept standard error
- Intercept_LCBxxPct
intercept xx% (default 95%) lower confidence bound
- Intercept_UCBxxPct
intercept xx% (default 95%) upper confidence bound
- Intercept_p_Value
intercept p-value
- Var_SiteInt
variance of the site intercepts
- Var_SiteTrend
variance of the site trends
- Corr_SiteIntSlope
correlation of site intercepts and site trends
- Var_Year
year variance
- Var_Residual
residual variance
- AIC
generalized Akaike Information Criterion
Details
For the simple linear regression (SLR) model, a design-based estimate of the
category proportion (categorical variables) or the mean (continuous
variables) is calculated for each time period (year). Four choices of
variance estimator are available for calculating variance of the design-based
estimates: (1) the local mean estimator, (2) the simple random sampling
estimator, (3) the Horvitz-Thompson estimator, and (4) the Yates-Grundy
estimator. For the Horvitz-Thompson and Yates-Grundy estimators, there are
three choices for calculating joint inclusion probabilities: (1) the Overton
approximation, (2) the Hartley-Rao approximation, and (3) the Brewer
approximation. The lm
function in the stats package is used to fit a
linear model using a formula
argument that specifies the proportion or
mean estimates as the response variable and years as the regressor variable.
For fitting the SLR model, the yearID
variable from the dframe
argument is modified by subtracting the minimum value of years from all
values of the variable. Parameter estimates are extracted from the object
returned by the lm
function. For the weighted linear regression (WLR)
model, the process is the same as the SLR model except that the inverse of
the variances of the proportion or mean estimates is used as the
weights
argument in the call to the lm
function. For the LMM
option, the lmer
function in the lme4 package is used to fit a linear
mixed-effects model for trend across years. For both the GLMM and LMM
options, the default Piepho and Ogutu (PO) model includes fixed effects for
intercept and trend (slope) and random effects for intercept and trend for
individual sites, where the siteID
variable from the dframe
argument identifies sites. Correlation between the random effects for site
intercepts and site trends is included in the model. Finally, the PO model
contains random effects for year variance and residual variance. For the GLMM
and LMM options, arguments cat_rhs
and cont_rhs
, respectively,
can be used to specify the right hand side of the model formula. Internally,
a variable named Wyear
is created that is useful for specifying the
cat_rhs
and cont_rhs
arguments. The Wyear
variable is
created by subtracting the minimum value of the yearID
variable from
all values of the variable. If argument invprboot
is FALSE
,
parameter estimates are extracted from the object returned by the lmer
function. If argument invprboot
is TRUE
, the boot
function in the boot package is used to generate bootstrap replicates using a
function named bootfcn
as the statistic
argument passed to the
boot
function. For each bootstrap replicate, bootfcn
calls the
glmer
or lmer
function, as appropriate, using the specified
model. design weights identified by the weight
argument for
the trend_analysis
function are passed as the weights
argument
for the boot
function, which specifies importance weights. Using the
design weights as the weights
argument ensures that bootstrap
replicates are representative of the survey population. Parameter estimates
are calculated using the object returned by the boot
function.
Author(s)
Tom Kincaid Kincaid.Tom@epa.gov
See Also
change_analysis
for change analysis
Examples
# Example using a categorical variable with three resource classes and a
# continuous variable
mydframe <- data.frame(
siteID = rep(paste0("Site", 1:40), rep(5, 40)),
yearID = rep(seq(2000, 2020, by = 5), 40),
wgt = rep(runif(40, 10, 100), rep(5, 40)),
xcoord = rep(runif(40), rep(5, 40)),
ycoord = rep(runif(40), rep(5, 40)),
All_Sites = rep("All Sites", 200),
Region = sample(c("North", "South"), 200, replace = TRUE),
Resource_Class = sample(c("Good", "Fair", "Poor"), 200, replace = TRUE),
ContVar = rnorm(200, 10, 1)
)
myvars_cat <- c("Resource_Class")
myvars_cont <- c("ContVar")
mysubpops <- c("All_Sites", "Region")
trend_analysis(
dframe = mydframe,
vars_cat = myvars_cat,
vars_cont = myvars_cont,
subpops = mysubpops,
model_cat = "WLR",
model_cont = "SLR",
siteID = "siteID",
yearID = "yearID",
weight = "wgt",
xcoord = "xcoord",
ycoord = "ycoord"
)