Regional Tropical Aboveground Biomass Mapping with L-Band Repeat-Pass Interferometric Radar, Sparse Lidar, and Multiscale Superpixels

: We introduce a multiscale superpixel approach that leverages repeat-pass interferometric coherence and sparse AGB estimates from a simulated spaceborne lidar in order to extend and from 57.1 Mg/ha to 46 Mg/ha. This work illustrates one of the numerous synergistic relationships between the spaceborne lidars, such as GEDI, with L-band SAR, such as PALSAR-2 and NISAR, in order to produce robust regional AGB in high biomass tropical regions.


Introduction
The sequestration of atmospheric carbon storage as Above Ground Biomass (AGB) is an important mitigating factor of the warming climate due to anthropogenic CO 2 [1,2]. Accurate, large-scale remote sensing of forests to monitor our global AGB stock is required for climate models and carbon policy [2,3]. L-band SAR is a powerful tool for monitoring Earth's AGB, particularly in lower biomass ecosystems, such as woodlands [3][4][5]. Techniques that are based on L-band backscatter alone suffer signal saturation at approximately 100-150 Mg/ha [6][7][8] and are impeded by environmental factors affecting vegetation's dielectric properties [7,9]. While new Polarimetric InSAR techniques, such as Random Volume over Ground [10,11] and Random Motion over Ground [12,13], promise to improve the characterization of forest structure and AGB in the future, there is currently no spaceborne instrument that is capable of providing required non-zero interferometric baseline measurements for these measurements. In anticipation of the upcoming NASA-ISRO SAR Mission (NISAR) mission, which is planned for launch in 2022, we seek to extend the capabilities of L-band products to estimate AGB across a larger AGB range, in particular dense tropical environments [14]. The NISAR science requirement is to estimate AGB within 20 Mg/ha RMSE with a spatial resolution of 1 ha for forest below 100 Mg/ha [6]. Although NISAR is an interferometric mission, albeit with zero-baseline repeat-pass orbit, current NISAR algorithms are based solely on the use of radar backscatter [6]. Here, we propose and evaluate the use of repeat-pass interferometric coherence to improve the estimation of AGB and extend the range of applicability beyond 100 Mg/ha. NISAR will operate in repeat-pass interferometric mode throughout its lifetime, monitoring surface deformation that results from earthquakes, volcanic activity, and other geologic and glacial processes. NISAR will repeat its orbit within 350 m to observe changes in phase that result from surface displacement rather than elevation or volumes. Nonetheless, this so-called zero spatial baseline interferometric measurement is related to forest structure [15][16][17], due to the temporal decorrelation of the repeat-pass interferometric signal [12,13,18]. Generally speaking, as noted in [15], the repeat-pass coherence is dominated by temporal decorrelation for shorter temporal baselines, barring significant weather or phenological events. In such cases, we expect areas with larger, more complex canopy structures to have lower coherence precisely due to the inevitable motion of such canopy structures between image retrievals.
Spaceborne lidar data from Ice, Cloud, and land Elevation Satellite's (ICESat) Geoscience Laser Altimeter System (GLAS) has been shown to provide invaluable calibration data to estimate global and continental-scale vegetation structure, such as canopy height [19,20] and AGB [21,22]. While GLAS was not originally designed to study forests, the new spaceborne lidar system, NASA's Global Ecosystem Dynamics Investigation (GEDI) [23] was designed to estimate global AGB stock within ±51.6 degrees latitudes [23,24]. The global high resolution sparse lidar AGB estimates from these satellites represent an extraordinary data set to calibrate high resolution regional vegetation models based on imaging instruments, such as Landsat [25,26], TanDEM-X [27][28][29], or L-/P-band SAR data [5,30].
In this paper, we investigate a methodology that is based on an L-band repeat-pass interferometric radar operating at HH and HV polarizations (similar to NISAR) to map AGB in tropical forests. The proposed approach calibrates a set of radar-to-biomass models with simulated spaceborne lidar samples, and utilizes interferometric coherence and multiscale superpixel-derived features from the L-band imagery to enhance AGB estimates. We represent spaceborne lidar data by sampling from our reference LVIS-derived AGB map, so that our shots and tracks are spaced in a manner similar to GEDI [23]; a simulation of end-of-mission GEDI or ICESat-2 precisely is beyond the scope of this work. We plan to apply this methodology to the available spaceborne lidar data in future work.

Canopy Cover Map
We only evaluate our methodology on forested areas as determined with the Hansen et al. canopy cover product [31]. This product describes the percent canopy within a 30 m Landsat pixel in the year 2000. We adopt the definition that a forested pixel has at least 30% canopy cover, which is in the generally accepted range of 10% to 30% of canopy cover [32,33]. We select a strict condition of forest to reduce falsely identified forest pixels, as suggested in [33]. Such a forest criterion has also been used by Hansen et al. for tropical forest analyses [31,34]. Our reported AGB ranges and accuracy metrics only consider pixels within this canopy cover mask.

AGB Reference Map
We validate our biomass estimation framework with an LVIS-derived AGB reference map [35]. Lidar-derived AGB is obtained from models for tropical lidar biomass models on the LVIS data, as described in [35][36][37]. Using these lidar-derived maps, we train our model along the hypothetical spaceborne tracks across our study area, in which along and across track spacing is similar to GEDI. We describe the simplified sampling strategy precisely in Section 3.2. In future work, we plan to apply this approach directly to the GEDI end-of-mission catalog. For the analysis of potential biomass errors from GEDI-backscatter fusion, we refer to [37].
In Table 1, we record the descriptive statistics of the AGB at each of the sites and sensors, including the coefficient of variation. Note that the two different L-band sensors will have different coverage at a particular site. Thus, the intersection of the SAR data with the reference AGB map will differ for each sensor and so too will the observed AGB range. In Figure 1, the intersection of the AGB reference data and the extents of the corresponding SAR image products are shown at each of the three sites. We note that Lope has the highest mean AGB, but Mondah has the highest coefficient of variation, given its lower mean AGB but comparable variance. In the table, we also note the descriptive statistics of the areas with AGB above and below 100 Mg/ha. From the last column of Table 1, Ogooue and Lope are predominantly composed of high AGB, whereas Mondah has more equally distributed low and high AGB areas. Table 1. Below, descriptive statistics of the LVIS-derived reference AGB in the study areas with respect to sensor coverage area. The µ AGB , µ AGB < 100 , and µ AGB > 100 are the mean AGB, the mean AGB for areas with AGB below 100 Mg/ha, and the mean AGB for areas with AGB above 100 Mg/ha, respectively. Similarly defined are σ AGB , σ AGB < 100 , and σ AGB > 100 for standard deviations over our study areas. The last column records the ratio of high AGB (above 100 Mg/ha) area to the total forest in the study area.

SAR Products
For each of our three sites, we separately estimate AGB for each L-band sensor: once using UAVSAR data and once using PALSAR-2 data. To do this, we obtain a backscatter image and a coherence image, as described in the next sections. Figure 2 shows a flowchart briefly summarizing our methodology from input SAR products and sparse lidar training data to the final AGB estimates.

Backscatter
We apply radiometric terrain correction (RTC) from [38] on dual-polarization (HH and HV) γ 0 backscatter images from PALSAR-2 and UAVSAR in order to remove the effects of topography on SAR backscatter. For both sensors, we download imagery in radar coordinates, i.e. slant-range products. While UAVSAR collected fully polarimetric data, we restrict ourselves to HH and HV channels as these have the greatest availability across the existing PALSAR-2 catalog and will be globally available with NISAR [39]. For UAVSAR, we use the multilooked slant-range products from the UAVSAR data page [40] and the RTC methodology described in [38]. For PALSAR-2, we used Level 1.1 slant-range products and subsequently geocode the imagery with the Gamma software package [41] to obtain similarly corrected backscatter images from the PALSAR-2 products. The RTC process removes the effects of viewing and target geometry on the observed backscatter imagery.
We despeckle each polarization of the RTC images to reduce erroneous speckle correlations with biomass and obtain superpixel segments that represent landscape features rather than noise-related shapes. We apply total variation (TV) denoising [42] within the homomorphic framework described in [43]. Total variation denoising removes zero mean additive Gaussian noise, preserving edges and discouraging localized oscillations [42,44]. We follow the image transformations, as in [43], so that the multiplicative noise from the L-band imagery is removed with this filter. Before TV denoising, we also perform a floor ( · ) and ceiling ( · ) of the top and bottom percentile of backscatter, respectively, to address the small number of pixels that were not successfully corrected due to layover and shadow in the backscatter imagery. We select the weight for total variation denoising by estimating the noise variance according to [45] and then selecting the regularization parameter according to [44,46]. Because images are frequently too large to hold in memory, we perform this despeckling in a 1000 × 1000 moving window with 980 stride, linearly interpolating between the two windows. We did not visually observe artifacts from this windowed denoising. We report the important backscatter statistics over canopied areas for each sensor in Table 2. We note that Mondah has the lowest mean backscatter, which is expected because of its lower biomass range (Table 1). We construct an RGB image while using the despeckled HH and HV backscatter images. We use HH, HV, and HH/HV as the red, green, and blue channels, respectively, as described in [47]. We regularize the polarization ratio HH/(HV + ε) with ε = 1 × 10 −5 to prevent division by 0. We then scale each channel to the interval [0, 1] for visualization. We remark that the polarization ratio HH/ HV has the physical interpretation of the ratio between non-volumetric and volumetric scattering [48]. In Figure 3, we compare the RTC image with the final despeckled product that we use for the subsequent estimation. Indeed, the despeckling preserves important boundaries between vegetated areas and river systems while removing the small local oscillations noticeable in the inland forest.
We reproject the backscatter into the reference frame of our lidar-derived AGB map with 50 m resolution. Therefore, the spatial resolution of our SAR products are reduced from 10 m to 50 m.

Coherence
Coherence is the magnitude of the correlation between two retrieved complex SAR images [49]. For UAVSAR, the temporal baselines range from a few hours (Mondah, Lope) to as long as eight days (Ogooue). For PALSAR-2, the temporal baseline is 14 days, which is the shortest available; we found that longer temporal baselines were not as useful for this correlative analysis due to higher levels of temporal decorrelation. We found only one 14-day PALSAR-2 pair for each of the tropical sites in the PALSAR-2 repository [50]. To obtain coherence with PALSAR-2, we used the open-source ISCE2 software package using the default parameters associated with their stripmap application [51] and used the 30 m Digital Elevation Model (DEM) from [52]. The window, in which the coherence is computed, is automatically determined with ICSE2, so that the resolution of the coherence image approximately equals the resolution of the DEM used during processing [51]. For UAVSAR, we utilized the open-source Kapok software [11] and applied Kapok to the coregistered image stacks available at the UAVSAR data page [40]. We use the window size from [53], which provides approximately 10 m resolution coherence images. We remark that the slant-range products used to generate the coherence and the backscatter imagery are different and, therefore, have different availability on [40]. We made sure that all of the products used in our analyses were captured within a window of a few weeks. We also note that both our original SAR products (coherence and backscatter) have finer spatial resolution than our reference AGB map. We perform the aforementioned processing of each SAR product at the finer resolution to ensure that our input SAR products that are projected to the reference AGB frame have high fidelity to the higher spatial resolution SAR products.

Generating Features from Multiscale Superpixels
Superpixels provide a coarse spatial segmentation of a scene into contiguous areas of homogeneous image intensity, in this case backscatter. Superixel analysis, also called "object-based" analysis in remote sensing [54,55], is applied to object detection [56,57], image classification [58,59], and change detection [60,61].
Here, we use superpixel segmentations of increasing size to generate features for our regional estimation model. Our method allows for us to characterize features of the SAR images at multiple scales simultaneously without the need to specify a scale a priori (e.g., [59,60]). We determine our segments using the backscatter RGB image from Section 2.3.1 while using the graph-based approach from [62]. The graph-based approach translates the image pixels into an 8-connected grid with each pixel as its own segment. Subsequently, the algorithm iteratively merges segments with its neighbors if a similarity criterion between segments is met in addition to ensuring that all segments have at least a user-defined minimum size [59,62]. To generate our multiscale features, we experimented with SLIC [63], an alternative superpixel methodology whose resulting segments have more uniform size than those derived from the graph-based approach. However, our AGB accuracy was almost always higher when using the graph-based segments. The superpixel features that are used in our study are the mean and standard deviation of the pixel values of each SAR product (coherence and backscatter) within the segments. The mean coherence features are particularly useful to smooth large spatial variability in observed coherence. A schematic of the multiscale features are shown in Figure 4, indicating how a pixel is related to these segmentations of increasing sizes. This multiscale approach was first proposed for object detection in [64] and it is similar to the image pyramid [65] applied in deep learning frameworks, such as [66]. Here, our generation of superpixel features is similar in spirit to those texture features generated by a 3 × 3 moving window in [14]. In Figure 5, the mean backscatter from the HV backscatter and HH coherence image within these multiscale superpixels are shown on a subset over Mondah. We also use standard deviation within a segment to derive features, though they are not shown. For the three superpixel scales, we selected minimum sizes 5, 25, 50, and fixed the scale parameter at 0.5. The size parameters were selected to ensure that the segments approximately doubled in mean size with each additional segmentation scale. We note the minimum size parameter is correlated with the mean segment size when the superpixel algorithm is applied to a fixed image, but a general relation across all imagery is not possible. The segmentation parameters used above are applied to all SAR backscatter imagery in this study. We did not use additional segmentation scales as the additional features did not markedly improve performance in our numerical investigations. We emphasize that we observed improvements in our AGB estimates with all sensors while using these fixed multiscale segmentation parameters.

Training and Validation: Simulating Regional Calibration with GEDI
We train our model on a sparse subset of pixels from our reference AGB map to approximate the coverage afforded with end-of-mission GEDI data [26,29,30], albeit our approach is highly simplified and does not represent the precise of spatial sampling of any existing spaceborne lidar. In our approach, we add tracks 120 m eastward modulo the image extents until we return to or pass our first track. We select pixels along lines oriented 51.2 clockwise below the equator starting at the corner of the reference AGB map spacing the sampled shots 60 m apart. We did not include the orthogonal tracks that will be available with GEDI, but note that the more data for training that we used, the better our model performed. We erred on the side of using less training data, as we expect GEDI data to be obscured by clouds in Gabon and, more generally, in tropical forests. We did not notice significant changes when translating the initial track across the site for different training sets; we observed no more than than 0.5 Mg/ha variation in the RMSE for such experimental translations of the training samples.  Figures  (e-h), are the analogous images derived using the mean HH coherence within superpixels. The maps are with respect to the reference AGB frame. The white areas in the imagery indicate where no data was collected by the lidar.
There are several advantages to the proposed regional training scheme, as noted in [4,9]. First, the training set contains a large range of AGB represented in the scene, which can help our model to overcome the model saturation that is typically associated with such AGB mapping [67]. Second, the calibration scheme is based on remotely sensed data and, therefore, overcomes the burden of collecting an adequate volume of field measurements in an area of interest. Thirdly, this calibration over the area of interest mitigates phenological and environmental conditions that could potentially degrade an L-band biomass model, particularly when the model is calibrated with respect to field measurements over a different area [9]. Although environmental factors can still degrade our AGB estimates, the calibration is done with respect to the relevant SAR imagery and the model adjusts accordingly. To expand on this last point, Ogooue's UAVSAR coherence proved to be less informative of AGB due to apparently uneven distributed rain over large parts of inland forest that occurred between the two acquisitions, causing temporal decorrelation in addition to its longer temporal baseline. However, the model training proceeded identically.
We note that we train our model on the entire sparse set of our lidar tracks whether or not they sample forest or non-forest areas, as defined in Section 2.1. We provide as many training points as possible over the region of interest because our sites are predominantly vegetated. If multiple land cover types are included in the area of interest, certain non-forest areas may need to be excluded during training to ensure our model can capture accurate correlative relationships between the SAR products and AGB. Our validation set is filtered according to the canopy mask from Section 2.1 in order to ensure that our errors are reported over strictly forested areas and not savanna or sparse low vegetation.
Random forest (RF) models have been successfully used for AGB estimation [14,68,69] and they are used in this work to estimate AGB [70] from the L-band features. The RF models can model non-linear relationships relating AGB to backscatter and handle multiple correlated inputs [70], which is especially important for our multiscale approaches in which features from segments over the same pixel will be correlated. The most relevant features can be determined post-training using impurity-based analysis with the GINI coefficient [70]. Additionally, the parallel implementation of RF models is efficient [71] and permits us to more easily train models with different inputs for comparison.
Our RF models ensemble 1000 trees with out-of-box scoring and an L 2 loss [70,71]. We experimented with binning our training data to simulate uniform sampling across the AGB range, but this did not improve the RMSE for any site or with any sensor. In Figure 6, we show a sample of the training and validation set over a small area in Mondah along with the estimates over the same test area.
We generate multiple realizations of the training data corrupted with additive Gaussian noise in order to assess potential impacts of uncertainty in GEDI's AGB products: where y is a biomass value and N (0, σ) is normally distributed with mean 0 and standard deviation σ. We ensure that y noisy ≥ 0 setting a numeric floor at 0. Increasing σ on the training data amplifies the size of errors introduced. We consider 100 different realizations of the training data according to Equation (1) for some fixed σ and consider the mean and variance of these errors to measure the noise impact on the final RMSE. This particular noise model is a highly simplified version of the errors that will be present in GEDI AGB estimates or other spaceborne lidars, as discussed in [72,73], and such noise simulation is beyond the scope of this work.

Results with a Low Biomass Mask
For each site and sensor, we use the corresponding reference biomass map to partition our area into areas that are above and below 100 Mg/ha. We apply our methodology to each area separately, obtaining a high and a low AGB model with a spatial resolution of 50 m (i.e., 0.25 ha). This will help to determine whether the results meet the NISAR AGB mapping requirement of less than 20 Mg/ha RMSE over areas with less than 100 Mg/ha [6], albeit with finer resolution. The partition of our study areas also allows us to more directly compare the results in [9], which only considers low AGB areas. Table 3 shows the results of this application. When we consider both high and low biomass areas separately over each site, the total nRMSE (here, we normalize according to the mean of the reference AGB) is <51%. Furthermore, our AGB maps have negligible bias, less than 1 percent of the true mean AGB in all cases. For areas below 100 Mg/ha, the AGB map in Lope and Mondah have RMSE below 20 Mg/ha at 50 m resolution, meeting the NISAR requirement. However, in Ogooue, the RMSE for both sensors slightly fell short of the NISAR RMSE requirement by 3 Mg/ha with UAVSAR and by 6 Mg/ha with PALSAR-2. For UAVSAR, the coherence appears to be degraded by rain that occurred between the repeat-pass radar acquisitions, causing significant temporal decorrelation not associated with vegetation structure, which might explain this poor performance. For areas that are above 100 Mg/ha, our RMSE is approximately 10 Mg/ha lower than the standard deviation of AGB for all sites and sensors, indicating that our model outperforms the mean estimator in this AGB regime. In other words, the model partially overcomes the model saturation typically associated with backscatter-AGB models [6,9]. In fact, for UAVSAR, our model's error is approximately 20 Mg/ha less than the standard deviation in the high AGB areas for all sites. As expected, the UAVSAR-based estimates were better than those from PALSAR-2, since UAVSAR has a higher signal-to-noise ratio, finer spatial resolution, and it has shorter temporal baseline (Table 2) than PALSAR-2.
The RMSE across all sites is more than double the RMSE in the study by Bouvet [9] that used a similarly derived reference AGB map. However, in [9], their validation area is primarily savanna and sparse woodlands. In contrast, our low AGB areas consist of short and dwarf Mangroves in Mondah and boundaries of dense tropical forest in Lope and Ogooue. Additionally, Lope and Ogooue are quite topographically varied ( Table 2). The methodology performs best over Mondah, which has the flattest terrain and the smallest AGB range of the three sites that we considered here.

Training on the Entire AGB Range
We also apply our estimation methodology to the full area without partitioning high and low AGB areas prior to training. This represents a more realistic application of the proposed methodology, in which a reference map is not available. Unsurprisingly, the accuracy, particularly on the low biomass areas suffers: the nRMSE for all areas is now <60% and the bias <3 Mg/ha for all sites and sensors ( Figure 7). The plots over Lope and Ogooue in Figure 7 appear to indicate a saturation of the AGB model between 300 and 400 Mg/ha. Although our AGB estimates frequently overestimate low AGB values and underestimate high AGB values, particularly in Lope and Ogooue, the overall low bias of our results indicates that our method provides an accurate regional estimate of aggregate AGB stock.
Several factors can explain the poorer model performance over low AGB areas. Firstly, our RF is trained using the L 2 loss and attempts to reduce absolute errors across the AGB range. These absolute errors are larger over the high AGB areas than the low AGB areas. In fact, since many of the sites are predominantly high AGB, accuracy in high AGB areas is more important to the RF model than accuracy in low AGB areas. In Table 4, we observe that the RMSE over the low AGB areas might be greater than or equal to the RMSE over the high AGB areas. Secondly, RFs are ensembles of regression trees that typically cluster estimates around the modes of the target distribution; in our case, this leads to underestimating high AGB values and overestimating low AGB values. Such errors have been observed in other AGB studies employing RFs [14,68]. Including more high quality optical imagery, such as Sentinel-2, for example, would allow our RF to better differentiate high and low AGB areas during training.  From the error metrics that are reported in Figure 7, we note, as in the previous section, that the UAVSAR estimates have lower absolute errors and higher model R 2 than those from PALSAR-2. The RMSE at Mondah is lowest when compared to the other sites and we conjecture, due, in part, to Mondah's flat terrain, lower AGB range and higher coefficient of variation. We also discuss in Section 4.2.1 that the PALSAR-2 coherence at Mondah is highly correlated with AGB.

Coherence and AGB
We show that the use of InSAR coherence can significantly improve the AGB estimates over those models that do not use InSAR coherence as an input. Given that these improvements were observed across all sites and sensors to varying degrees, utilizing the coherence for AGB estimation represents a significant opportunity for improving the range of observable AGB.
In Table 4, we compare our model across each site and sensor. For the comparison, we use all available superpixel features to obtain the lowest possible RMSE for both models. We found the most significant improvement at Mondah, where the addition of coherence reduces RMSE by approximately 9 Mg/ha for both sensors. The reduction in RMSE over high AGB areas was even larger, with PALSAR-2, the RMSE was reduced by nearly 14 Mg/ha over high AGB areas. Table 4. Below, we use the model with all superpixel features and compare the difference when the model uses backscatter with ("w/") and without ("w/o") coherence. In Figure 8, we plot the PALSAR-2 HH coherence against AGB within superpixel segments of minimum size 25 pixels (6.25 ha) with at least 30 percent of pixels being canopy, as specified by our forest mask from Section 2.1. We report the Spearman ρ there (−0.87), because there is a monotonically decreasing nonlinear correlation. This basic correlative analysis between coherence and AGB echoes those results from [18,28]. Unfortunately, the clear correlative relationship between coherence and AGB observed at Mondah with PALSAR-2 was not present elsewhere (in many cases, the Spearman correlation coefficient is not even defined). The apparent correlation of coherence with AGB at Mondah might be related to its high coefficient of variation of AGB, its lower AGB range, the flat terrain there, or some combination. We reemphasize that, in all of our numerical experiments, the inclusion of coherence reduced the RMSE of our AGB estimates.
We empirically illustrate that the inclusion of InSAR coherence improves our estimates, even when including noise to the training data using the noise model from Equation (1), because the GEDI mission's AGB estimates will still contain errors. In the top row of Figure 9, we compare the RF models trained with and without InSAR coherence when corrupted using this noise model. We select 50 equally spaced values of σ in [0, 200]. For each fixed σ, we generate training data corrupted by the additive noise of Equation (1) 100 times and apply the same training and validation methodology for each set of noisy training data. The figure shows the mean RMSE and one standard deviation of these errors for these two scenarios. As expected, the RMSE increases when our AGB model is trained on noisier training data. Our assertion that coherence reduces the RMSE of the AGB estimates remains valid for the majority of these models.

Multiscale Superpixels
The multiscale superpixel-derived features improve our regional AGB biomass estimates. Indeed, our results significantly improve across all sites and sensors as we add more superpixel-derived features as inputs; for UAVSAR, the improvement in RMSE is consistently over 12 Mg/ha and in some cases over 20 Mg/ha. Table 5 shows the results for the comparison of segment-derived features.
We visually compare the UAVSAR estimates that were made from models using exclusively pixelwise features and models using the additional three superpixel-derived features in subsets over Mondah (Figure 10), Lope (Figure 11), and Ogooue ( Figure 12). Visually, the results are more accurate using all of the superpixel features than those using just pixelwise features. However, the final estimates from the models using the superpixel-derived features still appear to be smoother than the reference AGB due to the spatial averaging used to generate such features. As in Section 4.2.1, we inspect the mean RMSE from 100 models trained on data generated using the noise model from Equation (1) for fixed σ in [0, 200] to compare models trained with just pixelwise data with those trained with all of the superpixel-derived features. This comparison is found in the second row of Figure 9. As expected, the RMSE increases with σ for both cases, but the models utilizing all the superpixel features consistently outperform those models that use only pixelwise data across the noise scenarios considered.

Inspecting Model Importances
We can inspect the feature importance of the RF using the method from [70]. Because our approach is non-parametric and not physically derived, looking at feature importances provides insight into this "black-box" machine learning approach. Roughly speaking, the most important features are those that when removed from the model input result in the highest decrease in estimation accuracy evaluated with respect to our training data. The feature importances are encoded in a vector of length that is equal to the number of inputs whose entries sum to 1 and whose relative value represents the corresponding feature's importance. In Table 6, we show the feature importances of a model trained on areas having AGB above 100 Mg/ha while using all of the possible superpixel-derived features. We label superpixel scales as 1, 2, 3 according to the minimum size 5, 25, 50 pixels used in our segmentation. We highlight that HH coherence is the most important feature at Mondah when using PALSAR-2 data, which is where we noticed the high Spearman correlation between AGB and coherence in Figure 8. Indeed, a coherence product at a particular scale appears in the top three features for all sites and sensors, except Lope imaged with UAVSAR; Lope had the largest fraction of high biomass and the lowest coefficient of variation with respect to biomass of the three sites considered and these could contribute to the fact that coherence was less useful there. We see that the polarization ratio (HH/HV) also appears as an important feature for numerous sites. This feature roughly represents the ratio of non-volumetric to volumetric scattering in canopy areas [48] and it can be used to distinguish dense forest and sparser vegetation structures.

Conclusions and Future Work
We introduced a regional above ground biomass (AGB) estimation framework leveraging sparse lidar training data and utilizing InSAR coherence and multiscale superpixels. We evaluated this framework over dense tropical forests in Gabon in order to demonstrate the efficacy of this approach for high AGB mapping while using NISAR-like data.
If we are able to partition our study area into low and high AGB areas, our regional AGB estimation methodology was able to meet NISAR RMSE requirement over low AGB areas in Lope and Mondah and fell slightly short of the requirement in Ogooue. We note that all of our AGB maps are at finer spatial resolutions than NISAR specifies. Although we could have more easily met the NISAR requirements using 1 ha resolution, we did not to pursue this here, because our multiscale superpixel approach would likely be superfluous on lower spatial resolution products, and our training and validation sets would be significantly reduced in size. Our proposed approach was partially able to overcome model saturation in the high AGB areas, particularly at Mondah while using both PALSAR-2 and UAVSAR, even at this higher resolution.
When we did not partition our study area into low and high AGB areas, the nRMSE increased, and in such cases was bounded by 60%. Additionally, the AGB maps no longer met the NISAR requirements in low AGB areas. We highlight that our AGB estimates all had low bias over the tropical sites providing reasonable estimates of the aggregate regional stock. We also note that Ogooue and Lope are predominantly high AGB areas, so our methodology might be suitable for AGB mapping in tropical forests with high AGB. We conjecture that including additional high quality image products, such as Sentinel-2, may allow our proposed estimation framework to better distinguish between high and low AGB areas during model training.
We also showed how the addition of InSAR coherence and multiscale features can significantly improve our AGB estimates. With PALSAR-2 InSAR coherence, the RMSE was reduced by over 8 Mg/ha (over 8% increase in nRMSE) over Mondah with both sensors. We observed improvements for all of our sites and sensors though none as significantly. Because PALSAR-2 InSAR coherence resembles what will be globally available with NISAR, the integration of such products for regional AGB retrieval represents a promising opportunity and will be investigated further. In future work, we hope to apply our framework to larger areas while using the GEDI mission data to calibrate our regional model. Funding: This study was primarily funded through NASA's Carbon Monitoring System (CMS, grant 15-CMS15-0055). Marc Simard was also funded by NASA Terrestrial and Ecology Program and a NISAR Science Team grant.

Acknowledgments:
The research was carried out at the Jet Propulsion Laboratory, California Institute of Technology, under a contract with the National Aeronautics and Space Administration. c 2020. California Institute of Technology. Government sponsorship acknowledged. We are grateful to the Agence Spatiale Gabonais d'Etudes et d'Observations Spatiales (AGEOS) for coordinating and overseeing the AfriSAR campaign making these Lidar-SAR analyses possible. We are grateful for JAXA's permission to use of ALOS-2/PALSAR-2 imagery over the tropical sites for this research purpose. We appreciate the helpful conversations with Creck Dassi, Nathan Thomas, and Daniel Jensen.

Conflicts of Interest:
The authors declare no conflict of interest.

Abbreviations
The following abbreviations are used in this manuscript: