Mangrove above-ground biomass and production are related to forest age at Low Isles, Great Barrier Reef

Above-ground biomass contributes a large proportion of mangrove carbon stock; however, spatio-temporal dynamics of biomass are poorly understood in carbonate settings of the Southern Hemisphere. This influences capacity to accurately project the effects of accelerating sea-level rise on this important carbon store. Here, above-ground biomass and productivity dynamics were quantified across mangrove age zones dominated by Rhizophora stylosa, spanning a tidal gradient atop a reef platform at Low Isles, Great Barrier Reef, Australia. Above-ground biomass was extrapolated across the forest using field plot data, allometry, a canopy height model derived from remotely piloted aircraft (RPA) LiDAR and regression analyses. Above-ground biomass production was calculated as mean annual biomass increments and canopy production was determined using RPA-derived multispectral imagery and a Normalised Difference Vegetation Index. Mangrove above-ground biomass was estimated at 519.7 ± 3.11 t ha-1 and increased with age up to the oldest forest (812.0 ± 12.9 t ha-1), believed to be ~135 ± 40 years old. Above-ground biomass was explained by age and tidal position (r2 >0.8), with a positive association between the two predictor variables. Above-ground biomass production peaked at lowest intertidal positions in the youngest forest aged <11 years at 36.3 t ha-1 yr-1, steadying thereafter, with a mean of 12.5 ± 5.4 t h-1yr-1 across the island. Production in the canopy remained high until the oldest forest and was negatively associated with age and tidal position (r2 >0.9). Declining production in the older zones corresponded to forest ageing, tidal positions becoming suboptimal for growth and increased exposure to prevailing winds and cyclones. By developing relationships between above-ground biomass accumulation and age and tidal position, this study informs parameterisation of models of the response of biomass to sea-level rise but requires additional information about relationships between substrate evolution and forest development and age.

Methods
Collection and analysis of field data

A total of 18 field plots (10 × 10 m) were established in the *R. stylosa-*dominant forest to collect tree measurements. Plots were positioned to span mangrove age zones and varying forest structures, assessed using historic mapping and field observations, respectively. Plots 1A–D and 2A–D are positioned in the 1848–1928 forest and plot 3 are in the 1928–1945 forest but close to the margin of the 1945–1973 forest. Plots 4A–D are in the 1945–1973 forest but also overlap with the 1973–2012 forest. Plots 5A and 5B are in the 1945–1973 forest and plots 6, 7 and 8 are in the 1973–2012 forest.

To calculate field plot AGB, tree trunk circumference at breast height (at 1.3 m) was measured and converted to DBH (cm). If trees were multi-stemmed, each stem was treated as a separate tree. For leaning trees, DBH was measured at 1.3 m along the natural height of the trunk, or above the highest prop root when this arose from the stem at a height above 1.3 m and this approach is consistent with the development of allometric equations from DBH. Individual tree heights were recorded in seven field plots but due to time constraints, one tree height was measured in plots 2A–D, a maximum tree height was measured in plots 1A–D and tree height data was not collected for the remaining five plots. Tree density was recorded as the number of individuals included in AGB calculations (i.e., above 1.3 m) and trees below this height were also tallied to determine shrub density. A total of 825 trees were measured for DBH and a total of 641 trees were measured for tree height. AGB for each tree was determined using Equation 1, an allometric equation based on DBH:

AGB (kg) = 10^A + B*log10DBH               (Equation 1)

Where: 10 represents the antilog of a base 10 log from rearranging the equation from (log Biomass = A + B log DBH) (B. F. Clough & Scott, 1989). The equation coefficients were derived from the total AGB of both R. stylosa and Rhizophora apiculata in northeastern Queensland (B. F. Clough & Scott, 1989), where A is -0.9789 and B is 2.6848. This equation was developed from trees with a DBH range of 3–23 cm, and therefore there is some uncertainty in values outside this range.

Collection and processing of remote sensing data

RPA remote sensing was used in this study as high-density LiDAR and high-resolution optical imagery was required and there was no other airborne LiDAR data available. An RPA survey of the mangrove forest was conducted in July 2023 using a DJI Matrice 300 RTK, fitted with a Zenmuse L1 and optical sensor to collect LiDAR and imagery in the visible spectrum (red, green, blue), respectively. Multispectral imagery (five-band) was collected in June 2022 using a DJI Phantom 4 fitted with a MicaSense RedEdge-M™ sensor. RPA surveys were conducted at low tides where possible to optimise detection of the substrate surface using LiDAR and spectral bands, both of which are sensitive to water (Klemas, 2013). Given the density of mangrove at this site, RPA-based LiDAR sensing was identified as the optimal approach for developing digital elevation models (DEMs), digital surface models (DSMs) and canopy height models (CHMs) due its high spatial resolution and precision.

Prior to the RPA survey, seven ground control points (GCPs) were placed across a range of elevations on the reef platform and within the mangrove forest. Trimble Real Time Extended Global Positioning System (RTX-GPS) points were also collected for easily identifiable features in the imagery such as microatolls and small mangroves on the reef platform. A total of 12 GCPs were surveyed using RTX-GPS (vertical error of 0.197 cm) to correct the vertical and horizontal offset of rasters.

Initial processing of the LiDAR data and the 2023 RPA imagery was conducted in DJI Terra (v.3.9.4), generating a las file and orthomosaic (2.78 cm/pixel), respectively. A multispectral orthomosaic (8 cm/pixel) was created by stitching images for each of the five spectral bands (blue, green, red, red edge, near-infrared) in Agisoft Metashape Professional** (Metashape v.1.6.3). Further image analysis was performed using ArcGIS Pro (v.3.3.1; ESRI) and R Studio (v.4.3.3; R Core Team 2024).

LiDAR ground points were classified to generate a DEM and DSM at 0.25 m resolution. The visible spectrum orthomosaic, DEM and DSM were corrected for horizontal offset by georeferencing using the GCPs and a 1st order polynomial correction (RMS errors <3 cm). A reliable RTX-GPS point was not available for the DJI Matrice 300 RTK base station at the time of the survey, which created a systematic vertical offset of the digital models. To rectify this, the average vertical offset between the DEM and GCPs was subtracted from both the DEM and DSM. The GCPs used for vertical correction were filtered to >0.28 m AHD to ensure submerged features did not influence the correction. This threshold was selected upon inspecting the orthomosaic for submerged features, as LiDAR returns do not effectively penetrate water.

To assess mangrove canopy greenness, the red and near-infrared multispectral orthomosaics were used to calculate NDVI.

Mangrove structure was quantified using a CHM (pixel size: 0.25 m), by subtracting the DEM from the DSM. The CHM was clipped to the 2023 mangrove extent and was resampled to a pixel size of 10 m (100 m2) to align with the field plot resolution for AGB. Field notes and GPS points guided the location and construction of square polygons representing the biomass plots (10 × 10 m). Calibration of the CHM was performed by linear regression analysis between field plot mean tree height where all tree heights were measured (n = 7) and the corresponding CHM pixel value to improve accuracy of the models. The calibrated CHM is hereafter referred to as Cal-CHM.

AGB was extrapolated across the mangrove forest using the Cal-CHM and a logistic regression relating field plot mean tree height and total plot AGB (n = 7) (Equation 2). As mangrove tree height asymptotes at a maximum height (Suwa et al. 2008), a logistic regression was most appropriate to explain the data. In some cases, linear regressions exhibit greater capacity to describe relationships between mangrove AGB and height (Suwa et al. 2008). However, they cannot effectively describe the AGB addition of very large trees that exhibit an asymptote in growth with age.

Y = c / 1 + e-a * (x + b)              (Equation 2)

Where: a is the growth rate estimated at 0.323, b is the inflection point estimated at 11.78 m, and c is the asymptote estimated at 31.77 t per plot. The field plot mean tree height was calculated for the seven plots with tree height data for all trees and is represented by x.

The Cal-CHM, AGB, DEM and NDVI rasters were resampled to 10 m resolution, clipped to the 2023 mangrove extent and interior ponds were removed to reduce outliers resulting from poor accuracy of LiDAR and multispectral imagery in submerged features.

Table. RPA survey details

RPA-multispectral
RPA-LiDAR and optical

Date
June 2022
July 2023

Drone, sensor
Phantom 4, MicaSense RedEdge-M™ sensor
DJI Matrice 300 RTK, Zenmuse L1, optical sensor

Flight speed (m s^-1^)

5

Altitude (m)
120
60

Conditions
Fine-cloudy, winds <15 km hr^-1^, low tide
Fine-cloudy, winds <15 km hr^-1^, low tide

Point cloud details

566.9 pts m^-2^ (spacing 0.042), 611,195,800 points, (16.4% ground, 83.6% undefined)

Optical image details
Red, green, blue, near infrared, red edge, 505 images per band (8 cm/pixel)
Red, green, blue, 970 images, (2.78 cm/pixel)

Table. DEM and DSM processing parameters.

DEM processing parameters
DSM processing parameters

LAS filter
Ground
All points

Interpolation type
Binning
Binning

Cell Assignment
Average
Maximum

Void Fill Method
Natural Neighbour
Linear

Cell size (m)
0.25
0.25