Habitat quality and rarity

Where is the forest of a landscape in good shape, and where is it worn down by the crops, roads and towns around it? InVEST's habitat quality model answers with a map of degradation and a map of quality, built only from the land use map, maps of the threats and two small tables. In this tutorial you will build it in Terraflux on Ribeirão Claro: you will describe three threats, see how each one degrades the habitat around it, choose the one constant that sets the scale of the quality, give a park its protection, compare the scenario of 2040 with 2020, and measure how rare each class has become since 1985. At the end you will know what the maps mean, how to set the model's parameters, and how not to over-read them.

The idea#

Habitat loss is easy to map: the forest that is no longer there. Habitat degradation is harder: the forest that remains but is hunted from the nearby road, invaded from the nearby pasture, burnt from the nearby fields, noisier and drier at its edges. The habitat quality model uses land use as a proxy for both. It asks three questions of every cell:

1. Is it habitat, and how good? Each class has a habitat suitability \(H_j\) from 0 (not habitat) to 1 (the best habitat), for the species or the biodiversity in general you have in mind.
2. Which threats are near, and how strong are they here? Each threat (crops, roads, a town) has a map, a maximum distance of influence, a relative weight, and a way its impact fades with distance: linearly or exponentially.
3. How sensitive is this habitat to each threat, and how protected is the cell? Each class has a sensitivity to each threat (0 to 1), and each cell an access (0 to 1) that lowers every threat where legal, social or physical barriers keep it out, as in a strictly protected park.

The model is one of InVEST's, described in its user's guide (Natural Capital Alliance, 2026). Terrado et al. (2016) adapted it to river basins and their freshwater habitats and compared it with field data: the quality it gave tracked the observed richness of vascular plants (R² = 0.76) and, less closely, the diversity of aquatic macroinvertebrates. Sallustio et al. (2017) mapped habitat quality over all of Italy at 20 m, with threats ranked by experts, and compared it with the network of protected areas.

A worked cell: a forest cell (sensitivity 0.6 to crops, 0.8 to roads) lies 500 m from the nearest cropland and 300 m from the nearest road; the town is beyond its reach. With this tutorial's threats (crops linear, 800 m, weight 0.7; roads exponential, 1,500 m, weight 0.8; urban exponential, 2,500 m, weight 1; the weights sum to 2.5), crops add \(0.28 \times (1 - 500/800) \times 0.6 = 0.063\) and roads \(0.32 \times e^{-2.99 \times 300/1500} \times 0.8 = 0.32 \times 0.550 \times 0.8 = 0.141\): \(D = 0.204\). With \(k = 0.35\) its quality is \(1 - 0.204^{2.5}/(0.204^{2.5} + 0.35^{2.5}) = 0.79\); with \(k = 0.05\) it would be 0.03.

Enlarge

Rarity#

With a base map (a land use of the past, "when intensive landscape management was relatively rare" in InVEST's words), the model also measures the rarity of each class: \(R_j = 1 - N_j / (N_{j,\text{base}} + N_j)\), where \(N_j\) is the area of class \(j\) on the current (or future) map and \(N_{j,\text{base}}\) its area on the base map. A class whose area has not changed has rarity 0.5; one that has shrunk has more than 0.5 (it has become rarer); one that has grown has less; a class absent from the base map has 0.

When to use it#

• As a proxy for biodiversity beside the other services, when there are no good maps of species: the scenarios of a plan compared for carbon, water and habitat at once.
• To find the remnants most at risk: habitat of high suitability and high degradation.
• To evaluate protection: how much the access of a park lowers the degradation of the habitat inside it.
• To track change: the quality of 2020 against a scenario, or against the past with the base map.

What it is not: it is not a model of any species' distribution, abundance or persistence. Its numbers are relative, set by the tables and by \(k\), and comparable between scenarios run with the same parameters, not between studies. It knows nothing about the size of patches or the connections between them; for those, see Fragmentation, MSPA and habitat networks. For one species, a species distribution model is the better tool.

What you need#

• The current land use map (lulc_2020.tif), projected in metres, integer codes.
• A map of each threat on the same grid, 1 where the threat is and 0 elsewhere: crops_2020.tif, roads_2020.tif, urban_2020.tif here. A threat can be a land use class (crops, urban) or something the land use map does not show (roads, mines, a railway).

Enlarge

• The threats, as a list in the operator's settings (below) or as InVEST's threats table, a CSV with the columns threat, max_dist (metres), weight (0 to 1) and decay (linear or exponential):

threat,max_dist,weight,decay
crops,800,0.7,linear
roads,1500,0.8,exponential
urban,2500,1.0,exponential

• The sensitivity table, a CSV by class with the columns lucode, habitat (\(H_j\), 0 to 1) and one column per threat, named as the threat, with the sensitivity of the class to it. Here forest is the best habitat (1), pasture and plantations poor habitat (0.3 and 0.2), the rest none:

• Optional: access polygons with an access attribute (here a park, 0, and its buffer zone, 0.5), a future map and its threat maps, a base map.

The values are illustrative. In a real study they come from the literature on the threats in your region and from experts (Sallustio et al. (2017) had experts rank the threats of each land use for all of Italy); document where each number comes from.

Build the model#

Step 1: the threat maps as one stack#

The operator takes the threat maps as a stack: several maps on one grid, each layer named. The names are what connects a map to its threat and to its column of the sensitivity table, so they must be the same in all three places.

format = "terraflux/1"

[[nodes]]
id = "crops"
op = "raster.read"
params = { path = "crops_2020.tif" }

[[nodes]]
id = "roads"
op = "raster.read"
params = { path = "roads_2020.tif" }

[[nodes]]
id = "urban"
op = "raster.read"
params = { path = "urban_2020.tif" }

[[nodes]]
id = "threat_maps"
op = "stack.of"
inputs = { crops = "crops.raster", roads = "roads.raster", urban = "urban.raster" }
note = "Layers named as the threats."

Step 2: the model#

The setting threats takes the list of threats, each with its four values. In the model file it reads:

format = "terraflux/1"

[[nodes]]
id = "lulc"
op = "raster.read"
params = { path = "lulc_2020.tif" }

[[nodes]]
id = "sensitivity"
op = "table.read"
params = { path = "sensitivity.csv" }
note = "lucode, habitat, crops, roads, urban."

[[nodes]]
id = "crops"
op = "raster.read"
params = { path = "crops_2020.tif" }

[[nodes]]
id = "roads"
op = "raster.read"
params = { path = "roads_2020.tif" }

[[nodes]]
id = "urban"
op = "raster.read"
params = { path = "urban_2020.tif" }

[[nodes]]
id = "threat_maps"
op = "stack.of"
inputs = { crops = "crops.raster", roads = "roads.raster", urban = "urban.raster" }

[[nodes]]
id = "habitat"
op = "services.habitat_quality"
inputs = { lulc = "lulc.raster", sensitivity = "sensitivity.table", threat_maps = "threat_maps.stack" }
params = { half_saturation = 0.05, threats = [
  { threat = "crops", max_dist = 800, weight = 0.7, decay = "linear" },
  { threat = "roads", max_dist = 1500, weight = 0.8, decay = "exponential" },
  { threat = "urban", max_dist = 2500, weight = 1.0, decay = "exponential" },
] }

Instead of the list, threats can be the name of InVEST's threats table: threats = "threats.csv". It gives the same maps (checked for this tutorial). Like the path of raster.read, the file name is taken from the model's folder, so the model and its table can move together.

Run it (Ctrl+Enter on the node). Its outputs:

• degradation: the total degradation \(D\) of each cell;
• threats: a stack with each threat's share of it, its layers named as the threats;
• quality: the habitat quality \(Q\);
• with a future or a base (step 4): degradation_future, quality_future, rarity, rarity_future, rarity_table, and with the base's threat maps degradation_base and quality_base.

Step 3: protection as access#

The park in the north is strictly protected and its buffer zone partly. Read access.gpkg with vector.read and connect its vector to the input access. Cells outside the polygons get access 1 (no protection). The input also takes a raster of access from 0 to 1, if you have one.

format = "terraflux/1"

[[nodes]]
id = "lulc"
op = "raster.read"
params = { path = "lulc_2020.tif" }

[[nodes]]
id = "sensitivity"
op = "table.read"
params = { path = "sensitivity.csv" }

[[nodes]]
id = "crops"
op = "raster.read"
params = { path = "crops_2020.tif" }

[[nodes]]
id = "roads"
op = "raster.read"
params = { path = "roads_2020.tif" }

[[nodes]]
id = "urban"
op = "raster.read"
params = { path = "urban_2020.tif" }

[[nodes]]
id = "threat_maps"
op = "stack.of"
inputs = { crops = "crops.raster", roads = "roads.raster", urban = "urban.raster" }

[[nodes]]
id = "access"
op = "vector.read"
params = { path = "access.gpkg" }
note = "A park (access 0) and its buffer zone (0.5)."

[[nodes]]
id = "habitat"
op = "services.habitat_quality"
inputs = { lulc = "lulc.raster", sensitivity = "sensitivity.table", threat_maps = "threat_maps.stack", access = "access.vector" }
params = { half_saturation = 0.35, threats = [
  { threat = "crops", max_dist = 800, weight = 0.7, decay = "linear" },
  { threat = "roads", max_dist = 1500, weight = 0.8, decay = "exponential" },
  { threat = "urban", max_dist = 2500, weight = 1.0, decay = "exponential" },
] }

Step 4: the future and the past#

For the scenario of 2040, read lulc_2040.tif and the three threat maps of 2040 into a second stack, and connect them to future and future_threat_maps; for the rarity, read lulc_1985.tif and connect it to base. The future needs its threat maps (the node says so if they are missing); the base does not (without base_threat_maps it only gives the rarity).

format = "terraflux/1"

[[nodes]]
id = "lulc"
op = "raster.read"
params = { path = "lulc_2020.tif" }

[[nodes]]
id = "future"
op = "raster.read"
params = { path = "lulc_2040.tif" }

[[nodes]]
id = "base"
op = "raster.read"
params = { path = "lulc_1985.tif" }

[[nodes]]
id = "sensitivity"
op = "table.read"
params = { path = "sensitivity.csv" }

[[nodes]]
id = "threats_2020"
op = "stack.of"
note = "crops_2020, roads_2020, urban_2020"

[[nodes]]
id = "threats_2040"
op = "stack.of"
note = "crops_2040, roads_2040, urban_2040"

[[nodes]]
id = "habitat"
op = "services.habitat_quality"
inputs = { lulc = "lulc.raster", sensitivity = "sensitivity.table", threat_maps = "threats_2020.stack", future = "future.raster", future_threat_maps = "threats_2040.stack", base = "base.raster" }
params = { half_saturation = 0.35, threats = "threats.csv" }

Read the results#

Degradation, threat by threat#

Enlarge

The degradation map is often the more useful of the two for management: it says which threat to act on and where. Its values have no units; they matter relative to each other and to \(k\).

Quality#

Enlarge

With the default \(k = 0.05\), nearly every forest cell is "degraded": the map says little. The reason is in the curve: \(k\) is the degradation at which quality is halved, and in this landscape the degradation of most forest is several times 0.05. InVEST's guide advises setting \(k\) to half of the highest degradation in the landscape: run once, read the largest value of degradation, and run again. Here the highest degradation is 0.70, so \(k = 0.35\).

Enlarge

Enlarge

Protection#

No access

Park and buffer zone

Access does not move the threats: a park whose forest is surrounded by crops keeps its quality only because you told the model the threats do not get in. Give access values you can defend (enforcement, fences, terrain), and show the result with and without them.

The scenario of 2040#

Enlarge

Enlarge

Two numbers describe a scenario better than one: the area of habitat and its mean quality. The quality-weighted area (each cell's quality times its area, summed) combines them into one and is a common summary in InVEST studies.

Rarity#

Enlarge

Rarity is a property of classes, not of places: every cell of a class has the same value. It is useful as a flag (which habitats are disappearing) and as a weight beside quality, not as a map of where to act.

Settings that change the answer#

The half-saturation constant#

\(k\) does not change the ranking of cells (a more degraded cell always has a lower quality), but it changes every number you will report:

Enlarge

PlayHalf-saturation constant k

Choose \(k\) once, from the current map, and keep it for the future and the past: a \(k\) recomputed for each scenario would hide the change you want to see.

The threats#

The maximum distances and the weights set how far and how strongly each threat reaches. Doubling the reach of the town (to 5 km) would carry its degradation over most of the map, which is 5.4 km from south to north; making roads linear instead of exponential would degrade more forest at medium distances. These are the numbers to justify from the literature (how far does hunting go from a road in your region? how far do edge effects of crops reach?) and to test: run the model with each changed by, say, half, and see which conclusions survive.

Common pitfalls#

• Names that do not match. The threat names must be the same in the threats (list or table), the layers of the stack, and the columns of the sensitivity table. A missing column stops the model with "the sensitivity table lacks the column(s) ...", a missing layer with "no current map of the threat ...".
• Threat maps with intensities. A road map with 0.3 for dirt roads and 1 for highways does not weaken the dirt roads in InVEST 3.20: both count fully. Make two threats instead (highways, dirt_roads) with their own weights.
• Weights out of range. Weights are from 0 to 1 and relative to each other; the model divides each by their sum.
• A maximum distance shorter than a cell. The model works in cells; a threat of 20 m on a 30 m grid reaches only its own cells.
• The default k. 0.05 is rarely right for your landscape; check the highest degradation.
• Comparing studies. Quality depends on the tables and on \(k\): a quality of 0.6 in one study and 0.8 in another say nothing about which landscape is in better shape.

Going further#

• Annual water yield, the next part.
• Carbon storage and sequestration for the same scenario in tons of carbon.
• Crop pollination by wild bees, part 8: the forest as the home of the bees that pollinate the farms around it.
• Ecological intactness, another index of how intact a landscape is, from the distance to human pressures.
• Fragmentation and MSPA for the shape of the habitat that remains.
• The operator's reference: services.habitat_quality, stack.of.

References#

• Natural Capital Alliance. (2026). InVEST 3.20.2 User's Guide: Habitat quality. https://storage.googleapis.com/releases.naturalcapitalproject.org/invest-userguide/latest/en/habitat_quality.html
• Sallustio, L., De Toni, A., Strollo, A., Di Febbraro, M., Gissi, E., Casella, L., Geneletti, D., Munafò, M., Vizzarri, M., & Marchetti, M. (2017). Assessing habitat quality in relation to the spatial distribution of protected areas in Italy. Journal of Environmental Management, 201, 129–137. https://doi.org/10.1016/j.jenvman.2017.06.031
• Terrado, M., Sabater, S., Chaplin-Kramer, B., Mandle, L., Ziv, G., & Acuña, V. (2016). Model development for the assessment of terrestrial and aquatic habitat quality in conservation planning. Science of the Total Environment, 540, 63–70. https://doi.org/10.1016/j.scitotenv.2015.03.064