add_linear_constraints  R Documentation 
Add constraints to a conservation planning problem to ensure that all selected planning units meet certain criteria.
## S4 method for signature 'ConservationProblem,ANY,ANY,character'
add_linear_constraints(x, threshold, sense, data)
## S4 method for signature 'ConservationProblem,ANY,ANY,numeric'
add_linear_constraints(x, threshold, sense, data)
## S4 method for signature 'ConservationProblem,ANY,ANY,matrix'
add_linear_constraints(x, threshold, sense, data)
## S4 method for signature 'ConservationProblem,ANY,ANY,Matrix'
add_linear_constraints(x, threshold, sense, data)
## S4 method for signature 'ConservationProblem,ANY,ANY,Raster'
add_linear_constraints(x, threshold, sense, data)
## S4 method for signature 'ConservationProblem,ANY,ANY,SpatRaster'
add_linear_constraints(x, threshold, sense, data)
## S4 method for signature 'ConservationProblem,ANY,ANY,dgCMatrix'
add_linear_constraints(x, threshold, sense, data)
x 

threshold 

sense 

data 

This function adds general purpose constraints that can be used to
ensure that solutions meet certain criteria
(see Examples section below for details).
For example, these constraints can be used to add multiple budgets.
They can also be used to ensure that the total number of planning units
allocated to a certain administrative area (e.g., country) does not exceed
a certain threshold (e.g., 30% of its total area). Furthermore,
they can also be used to ensure that features have a minimal level
of representation (e.g., 30%) when using an objective
function that aims to enhance feature representation given a budget
(e.g., add_min_shortfall_objective()
).
An updated problem()
object with the constraints added to it.
The linear constraints are implemented using the following
equation.
Let I
denote the set of planning units
(indexed by i
), Z
the set of management zones (indexed by
z
), and X_{iz}
the decision variable for allocating
planning unit i
to zone z
(e.g., with binary
values indicating if each planning unit is allocated or not). Also, let
D_{iz}
denote the constraint data associated with
planning units i \in I
for zones z \in Z
(argument to data
, if supplied as a matrix
object),
\theta
denote the constraint sense
(argument to sense
, e.g., <=
), and t
denote the constraint
threshold (argument to threshold
).
\sum_{i}^{I} \sum_{z}^{Z} (D_{iz} \times X_{iz}) \space \theta \space t
The argument to data
can be specified using the following formats.
data
as character
vectorcontaining column name(s) that
contain penalty values for planning units. This format is only
compatible if the planning units in the argument to x
are a
sf::sf()
or data.frame
object. The column(s) must have numeric
values, and must not contain any missing (NA
) values.
For problems that contain a single zone, the argument to data
must
contain a single column name. Otherwise, for problems that
contain multiple zones, the argument to data
must
contain a column name for each zone.
data
as a numeric
vectorcontaining values for
planning units. These values must not contain any missing
(NA
) values. Note that this format is only available
for planning units that contain a single zone.
data
as a matrix
/Matrix
objectcontaining numeric
values
that specify data for each planning unit.
Each row corresponds to a planning unit, each column corresponds to a
zone, and each cell indicates the data for penalizing a planning unit
when it is allocated to a given zone.
data
as a terra::rast()
objectcontaining values for planning
units. This format is only
compatible if the planning units in the argument to x
are
sf::sf()
, or terra::rast()
objects.
If the planning unit data are a sf::sf()
object,
then the values are calculated by overlaying the
planning units with the argument to data
and calculating the sum of the
values associated with each planning unit.
If the planning unit data are a terra::rast()
object, then the values
are calculated by extracting the cell
values (note that the planning unit data and the argument to data
must
have exactly the same dimensionality, extent, and missingness).
For problems involving multiple zones, the argument to data
must
contain a layer for each zone.
See constraints for an overview of all functions for adding constraints.
Other constraints:
add_contiguity_constraints()
,
add_feature_contiguity_constraints()
,
add_locked_in_constraints()
,
add_locked_out_constraints()
,
add_mandatory_allocation_constraints()
,
add_manual_bounded_constraints()
,
add_manual_locked_constraints()
,
add_neighbor_constraints()
## Not run:
# load data
sim_pu_raster < get_sim_pu_raster()
sim_features < get_sim_features()
# create a baseline problem with minimum shortfall objective
p0 <
problem(sim_pu_raster, sim_features) %>%
add_min_shortfall_objective(1800) %>%
add_relative_targets(0.2) %>%
add_binary_decisions() %>%
add_default_solver(verbose = FALSE)
# solve problem
s0 < solve(p0)
# plot solution
plot(s0, main = "solution", axes = FALSE)
# now let's create some modified versions of this baseline problem by
# adding additional criteria using linear constraints
# first, let's create a modified version of p0 that contains
# an additional budget for a secondary cost dataset
# create a secondary cost dataset by simulating values
sim_pu_raster2 < simulate_cost(sim_pu_raster)
# plot the primary cost dataset (sim_pu_raster) and
# the secondary cost dataset (sim_pu_raster2)
plot(
c(sim_pu_raster, sim_pu_raster2),
main = c("sim_pu_raster", "sim_pu_raster2"),
axes = FALSE
)
# create a modified version of p0 with linear constraints that
# specify that the planning units in the solution must not have
# values in sim_pu_raster2 that sum to a total greater than 500
p1 <
p0 %>%
add_linear_constraints(
threshold = 500, sense = "<=", data = sim_pu_raster2
)
# solve problem
s1 < solve(p1)
# plot solutions s1 and s2 to compare them
plot(c(s0, s1), main = c("s0", "s1"), axes = FALSE)
# second, let's create a modified version of p0 that contains
# additional constraints to ensure that each feature definitely has
# at least 8% of its overall distribution represented by the solution
# (in addition to the 20% targets which specify how much we would
# ideally want to conserve for each feature)
# to achieve this, we need to calculate the total amount of each feature
# within the planning units so we can, in turn, set the constraint thresholds
feat_abund < feature_abundances(p0)$absolute_abundance
# create a modified version of p0 with additional constraints for each
# feature to specify that the planning units in the solution must
# secure at least 8% of the total abundance for each feature
p2 < p0
for (i in seq_len(terra::nlyr(sim_features))) {
p2 <
p2 %>%
add_linear_constraints(
threshold = feat_abund[i] * 0.08,
sense = ">=",
data = sim_features[[i]]
)
}
# overall, p2 could be described as an optimization problem
# that maximizes feature representation as much as possible
# towards securing 20% of the total amount of each feature,
# whilst ensuring that (i) the total cost of the solution does
# not exceed 1800 (per cost values in sim_pu_raster) and (ii)
# the solution secures at least 8% of the total amount of each feature
# (if 20% is not possible due to the budget)
# solve problem
s2 < solve(p2)
# plot solutions s0 and s2 to compare them
plot(c(s0, s2), main = c("s1", "s2"), axes = FALSE)
# third, let's create a modified version of p0 that contains
# additional constraints to ensure that the solution equitably
# distributes conservation effort across different administrative areas
# (e.g., countries) within the study region
# to begin with, we will simulate a dataset describing the spatial extent of
# four different administrative areas across the study region
sim_admin < sim_pu_raster
sim_admin < terra::aggregate(sim_admin, fact = 5)
sim_admin < terra::setValues(sim_admin, seq_len(terra::ncell(sim_admin)))
sim_admin < terra::resample(sim_admin, sim_pu_raster, method = "near")
sim_admin < terra::mask(sim_admin, sim_pu_raster)
# plot administrative areas layer,
# we can see that the administrative areas subdivide
# the study region into four quadrants, and the sim_admin object is a
# SpatRaster with integer values denoting ids for the administrative areas
plot(sim_admin, axes = FALSE)
# next we will convert the sim_admin SpatRaster object into a SpatRaster
# object (with a layer for each administrative area) indicating which
# planning units belong to each administrative area using binary
# (presence/absence) values
sim_admin2 < binary_stack(sim_admin)
# plot binary stack of administrative areas
plot(sim_admin2, axes = FALSE)
# we will now calculate the total amount of planning units associated
# with each administrative area, so that we can set the constraint threshold
# since we are using raster data, we won't bother explicitly
# accounting for the total area of each planning unit (because all
# planning units have the same area in raster formats)  but if we were
# using vector data then we would need to account for the area of each unit
admin_total < Matrix::rowSums(rij_matrix(sim_pu_raster, sim_admin2))
# create a modified version of p0 with additional constraints for each
# administrative area to specify that the planning units in the solution must
# not encompass more than 10% of the total extent of the administrative
# area
p3 < p0
for (i in seq_len(terra::nlyr(sim_admin2))) {
p3 <
p3 %>%
add_linear_constraints(
threshold = admin_total[i] * 0.1,
sense = "<=",
data = sim_admin2[[i]]
)
}
# solve problem
s3 < solve(p3)
# plot solutions s0 and s3 to compare them
plot(c(s0, s3), main = c("s0", "s3"), axes = FALSE)
## End(Not run)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.