Simulation-Based Evaluation of the Estimation Methods of Far-Red Solar-Induced Chlorophyll Fluorescence Escape Probability in Discontinuous Forest Canopies

Abstract

The escape probability of Solar-induced chlorophyll fluorescence (SIF) can be remotely estimated using reflectance measurements based on spectral invariants theory. This can then be used to correct the effects of canopy structure on canopy-leaving SIF. However, the feasibility of these estimation methods is untested in heterogeneous vegetation such as the discontinuous forest canopy layer under evaluation here. In this study, the Discrete Anisotropic Radiative Transfer (DART) model is used to simulate canopy-leaving SIF, canopy total emitted SIF, canopy interceptance, and the fraction of absorbed photosynthetically active radiation (fAPAR) in order to evaluate the estimation methods of SIF escape probability in discontinuous forest canopies. Our simulation results show that the normalized difference vegetation index (NDVI) can be used to partly eliminate the effects of background reflectance on the estimation of SIF escape probability in most cases, but fails to produce accurate estimations if the background is partly or totally covered by vegetation. We also found that SIF escape probabilities estimated at a high solar zenith angle have better estimation accuracy than those estimated at a lower solar zenith angle. Our results show that additional errors will be introduced to the estimation of SIF escape probability with the use of satellite products, especially when the product of leaf area index (LAI) and clumping index (CI) was underestimated. In other results, fAPAR has comparable estimation accuracy of SIF escape probability when compared to canopy interceptance. Additionally, fAPAR for the entire canopy has better estimation accuracy of SIF escape probability than fPAR for leaf only in sparse forest canopies. These results help us to better understand the current estimation results of SIF escape probability based on spectral invariants theory, and to improve its estimation accuracy in discontinuous forest canopies.

1. Introduction

Solar-induced chlorophyll fluorescence (SIF) is emitted from within the photosynthetic apparatus of all higher plants and is a good indicator of carbon assimilation and plant physiological status [1,2]. Hence, SIF has been used to characterize the spatiotemporal dynamics of gross primary productivity (GPP), and strong correlations have been found between satellite-derived SIF and landscape-scale GPP across different biomes [3,4,5,6,7,8].

Many research studies have reported a strong linear relationship between TOC (Top of Canopy) SIF and ecosystem-scale GPP over different biomes [9]. However, whether or not a universal linear relationship exists across biomes at ecosystem-scale is still an open question [6,10]. Usually, GPP is calculated as the difference between NEE (Net Ecosystem Exchange) and ER (Ecosystem Respiration) [11], while the satellite-observed TOC SIF is calculated based on the directional canopy spectrum [12,13,14]. TOC SIF is only a part of the total emitted SIF of vegetation and is influenced by canopy structure and solar-viewing geometry. These factors weaken the relationship between SIF and GPP. For this reason, the total emitted SIF of vegetation (SIF emission from all leaves within the canopy) is more directly related to GPP than TOC SIF [13,14].

Canopy radiative transfer models (RTMs) offer an easy and direct way to estimate the total emitted SIF through model inversion. This process amounts to the compensation of the observed TOC SIF for the structure-dependent multiple scattering and re-absorption events that occur from the leaves to sensor. Today, several SIF models have already been developed that are capable of this, including the one-dimensional (1-D) SIF RTM, such as the Soil-Canopy Observation of Photochemistry and Energy (SCOPE) balance model [15] and three-dimensional (3-D) SIF RTM, such as the FluorWPS model [16], the FluorFLIGHT model [17], and the Discrete Anisotropic Radiative Transfer (DART) model [18]. However, a lot of input parameters from leaf to canopy are needed for these models to carry out the inversion process, and these input parameters are not always available on a global scale.

Spectral invariants theory (SIT) was developed and introduced to remote sensing to simplify the radiative transfer processes among leaf, shoot, and canopy [19,20,21]. When compared with RTMs, the radiation scattering processes within canopies can be characterized just by leaf albedo and three spectral invariant parameters: canopy interceptance (i_o), directional escape probability, and recollision probability under the assumption that soil background has no (or negligible) effects on canopy reflectance. Recently, a semi-empirical approach has been developed to estimate SIF escape probability by Yang et al. [22]. Here canopy SIF scattering at far-red band is estimated using TOC pure vegetation reflectance, in addition to two other parameters, canopy interceptance and leaf albedo. These parameters are combined to estimate the SIF escape probability based on SIT ($f_{yang}^{esc}=NI{R}_v/i_o·{\omega }_{\lambda }$, $NI{R}_v$ is defined as the vegetation canopy reflectance by vegetation, $i_o$ is the canopy interceptance, and ${\omega }_{\lambda }$ is the leaf albedo at the specific wavelength), which is finally used to derive total emitted SIF from observed TOC SIF.

As stated above, canopy interceptance (i_o) is required for the estimation of SIF escape probability and not commonly available on a global scale. To address this problem, Liu et al. [23] proposed a method to estimate the SIF escape probability using random forest with only TOC SIF and reflectance data. As this estimation method is based on an empirical approach, a more generalizable physical based approach is also sought. Zeng et al. [12] proposed an extended method for the estimation of SIF escape probability based on the concept of Yang et al. [22] (i.e., $f_{yang}^{esc}=NI{R}_v/i_o·{\omega }_{\lambda }$), fraction of absorbed photosynthetically active radiation (fAPAR) was used to replace the canopy interceptance for remote sensing applications in this method. In addition, the normalized difference vegetation index (NDVI) was used in this study to assess the reflectance of pure vegetation (NIR_v ≈ NIR_T$·$NDVI) to solve the “black soil problem” [24]. Zhang et al. [14] proposed to directly derive the canopy interceptance from remotely sensed data at hand (leaf area index (LAI) and clumping index (CI) satellite products) using the method of [25] (i.e., $i_{\circ }=1-\exp \left(-G\left(\theta \right)·LAI·CI/COS\left(\theta \right)\right)$).

However, these estimation methods were developed for homogeneous vegetation which tend to the theoretical and have been mainly validated using 1-D radiative transfer models. Their performance for the estimation of SIF escape probability has not yet been tested and evaluated for sparse vegetation canopies which typify large areas of the boreal forest [26], most such forests are discontinuous in structure and contain overstory and understory layers with varying species and phenologies which depend on factors such as site fertility [27].

In this study, we investigate whether the estimation methods of far-red SIF escape probability based on SIT holds for discontinuous forest canopies using the model-based assessment strategy. To achieve this, we used a 3-D RT model, DART, to simulate the TOC SIF and per-leaf emitted SIF (i.e., total emitted SIF) based on the 3-D radiative budget simulation results. After that, the escape probability of SIF in the forest area was directly calculated based on the 3-D simulation results. Subsequently, the performance of estimation methods for escape probability based on SIT was analyzed under different vegetation conditions, such as different leaf optical properties, canopy structure, and background optical properties.

2. Materials and Methods

2.1. Spectral Invariants Theory and SIF Escape Probability

The original estimation method (Equation (1)) for SIF escape probability at far-red band based on SIT was developed by [22]. It combines top of canopy reflectance NIR_v, canopy interceptance i_o, and leaf single scattering albedo ω_λ as

$$f^{esc}=\frac{NI{R}_v}{i_{°}·{\omega }_{\lambda }}.$$

(1)

In a remote sensing application, NIR_v is generally replaced by NIR_T or NIR_T$·$NDVI, where NDVI is used to solve the “black soil problem” [12]. In addition, fAPAR can be used to replace canopy interceptance i_o because fAPAR is simpler to estimate. Based on these considerations, we decided to calculate and evaluate SIF escape probability ($f_{SIT}^{esc}$: f^esc based on SIT) using four formulas shown in

Table 1. Four estimation formulas for Solar-induced chlorophyll fluorescence (SIF) escape probabilities based on spectral invariants theory (SIT) ($f_{SIT}^{esc}$).

$${\mathit{f}}_{\mathit{S}\mathit{I}\mathit{T}}^{\mathit{e}\mathit{s}\mathit{c}}$$	NIRT	$${\mathit{NIR}}_{\mathit{v}}\approx {\mathit{NIR}}_{\mathit{T}}·\mathbf{NDVI}$$
io	$$\frac{NI{R}_T}{i_{°}·{\omega }_{\lambda }}$$	$$\frac{NI{R}_v}{i_{°}·{\omega }_{\lambda }}$$
fAPAR	$$\frac{NI{R}_T}{\mathrm{fAPAR}·{\omega }_{\lambda }}$$	$$\frac{NI{R}_V}{\mathrm{fAPAR}·{\omega }_{\lambda }}$$

.

2.2. DART Simulation

2.2.1. Simulation of SIF Escape Probability

The Discrete Anisotropic Radiative Transfer (DART) model is a 3-D RT model [28] that works with realistic representations of plant architecture described by the foliage distribution, including leaf clumping at branch and crown levels, along with detailed geometry of stems and branches. Recently, the leaf fluorescence FLUSPECT model has been integrated into the DART model to simulate fluorescence from leaf to canopy level [18]. Besides the simulation of spectral total and SIF radiance that escape the canopy, DART also simulates the 3-D, 2-D, and 1-D spectral radiative budgets of vegetation scenes per type of scene element (e.g., ground, woody elements, and leaves). It can also simulate the radiative budget of individual facets (i.e., triangles) that are used to represent scene elements such as leaves. These radiative budgets are computed at different stages of the radiative transfer modeling: after direct sun illumination, after atmosphere illumination, and after each subsequent DART iteration. They include radiation that is intercepted, scattered, absorbed, and emitted (i.e., SIF in the short waves and thermal emission in the long waves), and also the irradiance and exitance at each face of each voxel with which the landscape is simulated. Units can be W/m²/µm or percentage of top of canopy irradiance. Therefore, DART can compute the light regime and instantaneous incoming photosynthetically active radiation (PAR) and fAPAR of each leaf element and each voxel and layer of the canopy.

The spectral SIF escape probability (f^esc) was defined as the ratio of canopy-leaving SIF observed by sensors (SIF_obs) to total emitted SIF from whole-canopy (SIF_tot) at a given wavelength [12],

$$f_{DART}^{esc}=\frac{SI{F}_{obs}}{SI{F}_{tot}},$$

(2)

where SIF_obs is the TOC SIF exitance. It corresponds to the angular integration of DART TOC SIF radiance images; more details about the simulation process can be found in [29]. SIF_tot is the sum of all SIF photons emitted by all leaves within the canopy in all directions [12,24]. SIF_tot was derived from the DART 3-D radiative budgets of SIF emitted in all directions for each leaf within the 3-D vegetation canopy.

2.2.2. Simulation of Canopy Interceptance

Canopy interceptance (i_o) is defined as the probability that an incoming solar photon interacts in first-order with leaves (or needles) within the canopy at any wavelength [21]. Because it is equal to one minus directional gap fraction, it is often derived from the simulation of directional gap fraction as in the works of [12,22] using the SCOPE model. In DART, i_o is the “Intercepted” term of the 3-D and 1-D vegetation radiative budgets at the end of the direct sun illumination stage (i.e., “Illudir” stage in DART), before the scattering stage. It is equal to 1—“Intercepted” terms of the ground and trunk/branch radiative budgets of the Illudir stage.

$$i_o^{DART}=\frac{R{B}_{leaf}^{intercepted}\left(\lambda \right)}{R{B}_{leaf}^{intercepted}\left(\lambda \right)+R{B}_{ground}^{intercepted}\left(\lambda \right)},$$

(3)

where $R{B}_{leaf}^{intercepted}$ and $R{B}_{ground}^{intercepted}$ are the leaf absorbed irradiance and ground intercepted radiation at a specific wavelength, respectively.

2.2.3. Simulation of fAPAR

fAPAR_canopy was defined as the fraction of absorbed photosynthetically active radiation (APAR) by the entire canopy (i.e., leaves and wood elements without ground) over the PAR region from 400 nm to 700 nm [30]. Accordingly, fAPAR_leaf was defined as the fraction of absorbed photosynthetically active radiation (APAR) only by leaves. To evaluate the impact of woody elements on the simulation and estimation of SIF escape probabilities, two types of fAPAR (i.e., fAPAR_canopy = fAPAR_leaf+woody and fAPAR_leaf) were used to estimate SIF escape probabilities (Section 3.4).

In DART, fAPAR_canopy and fAPAR_leaf were derived from the radiative budget of facets per type of scene element (i.e., leaf, wood element, ground).

$$\mathrm{fAPAR}=\frac{{{\displaystyle \int }}_{400 nm}^{700 nm}R{B}^{absorbed}\left(\lambda \right)d\lambda }{{{\displaystyle \int }}_{400 nm}^{700 nm}E\left(\lambda \right)d\lambda },$$

(4)

where $E\left(\lambda \right)$ is the incident irradiance, and $\lambda$ is wavelength over the spectral region from 400 to 700 nm, $R{B}^{absorbed}\left(\lambda \right)$ is the irradiance absorbed by leaf only (i.e., $R{B}_{leaf}^{absorbed}\left(\lambda \right)$) or entire canopy (i.e., $R{B}_{canopy}^{absorbed}\left(\lambda \right)$ = $R{B}_{leaf}^{absorbed}\left(\lambda \right)$ + $R{B}_{woody}^{absorbed}\left(\lambda \right)$). In addition, $R{B}^{absorbed}\left(\lambda \right)$ was calculated by area-weighted average as below:

$$R{B}^{absorbed}\left(\lambda \right)={{\displaystyle \sum }}_{i=1}^n\left(R{B}_i^{absorbed}·S_i^{element}\right)/S_{canopy}^{element},$$

(5)

where n is the number of leaves in the canopy, $R{B}_i^{absorbed}$ and $S_i^{element}$ are the absorbed irradiance and area of scene elements of type i (e.g., leaf, woody element), respectively, and $S_{canopy}^{element}$ is the sum of areas of all elements within the canopy.

2.2.4. Simulation Experiment

We conducted a numerical simulation experiment to simulate SIF escape probability using the 3-D DART model and then using it to evaluate the validity of the remote estimation methods described in Section 2.1. DART provides a realistic representation of plant architecture [18,29]. Hence, DART provides a way to simulate and evaluate SIF escape probability in a more realistic vegetation scene relative to simpler models, especially for sparse forest stands which are found across the boreal forest [31].

A 3-D Silver birch tree (Betula pendula Roth) reconstructed from Terrestrial laser scanning (TLS) data, measured in the vicinity of the Station for Measuring Ecosystem-Atmosphere Relations II (SMEAR II) in southern Finland (24.31478° E, 61.84335° N), was used to generate a 3-D simulation scene. We used an option of repetitive scene of the DART model to eliminate boundary effects [29]. The scene horizontal dimension was set to 4 m for default configuration (i.e., fractional vegetation cover (FVC) = 0.6, see

Table 2. Main input parameters of the Discrete Anisotropic Radiative Transfer (DART) simulations with indication of the default configuration.

	Parameters	Unit	Values	Default
Leaf optical properties	Chlorophyll content (Cab)	μg/cm2	15, 35, 60	35
Canopy structure	Fractional vegetation cover (FVC)	[-]	0.3, 0.6, 0.8	0.6
Background spectral properties	Background reflectance spectra	[-]	Three reflectance spectra (see Figure S1: Ref-1, 2, and 3) and a non-reflecting background	Ref-2
Solar-view geometry	Solar zenith angle (SZA)	[°]	20–70° with a step of 10°	30

). By adjusting the scene horizontal dimensions, the fractional vegetation cover (FVC) was modified from sparse forest (FVC = 0.3) to a dense forest (FVC = 0.8).

Due to the computational constraints, four key input parameters of the DART model were varied to generate the simulation database (see Table 2), including leaf chlorophyll content (Cab), fractional vegetation cover (FVC), background spectra, and solar zenith angle (SZA). Variable values of leaf chlorophyll content were used to represent the variation of leaf optical properties. Varying values of fractional vegetation cover were used to describe the canopy structure of different forest types: 0.3, 0.6, and 0.8 are referred to as very sparse, little sparse, and dense forest, respectively. Given the two-layer structure of the forest and its seasonal variation in spectral properties, three typical background spectra were adapted to represent the change of background spectral properties. In addition, a non-reflecting background was used to test the ability of SIF escape probability that was estimated using SIT. Contrary to the DART simulations carried out by previous studies in which the solar zenith angle was set to a constant value [12,23], the solar zenith angle was varied here from 20° to 70° with a step of 10° to represent the seasonal variation for most of the boreal forest areas. Details about the other default variables for the DART simulations can be found in [29].

All simulations were performed between 640 and 850 nm at a spectral resolution of 7.5 nm. Directional and total emitted SIF at 760 nm were used for the simulation and estimation of SIF escape probability. All canopy reflectance and SIF were retrieved in nadir viewing direction (i.e., view zenith angle equal to 0) to keep consistent with the default observation mode of current SIF satellite (such as the OCO-2 satellite). NDVI was derived from the canopy reflectance in the red (650 nm) and near-infra-red (NIR) (800 nm) bands according to [32]. SIF escape probabilities simulated by the DART model ($f_{DART}^{esc}$) were regarded as ‘real’ value, the relative estimation accuracy of SIF escape probabilities using SIT ($f_{SIT}^{esc}$) were calculated as ($f_{SIT}^{esc}-f_{DART}^{esc}$)/$f_{DART}^{esc}$.

2.3. Remote Estimation of Canopy Interceptance

In addition to canopy reflectance, canopy interceptance (i_o) is a crucial parameter for the estimation of SIF escape probability [22]. Although i_o is not easily accessed directly from remote sensing measurements, it can be indirectly estimated from remotely sensed data using equations such that of [25]:

$$i_{\circ }=1-\exp \left(-G\left(\theta \right)·LAI·CI/COS\left(\theta \right)\right),$$

(6)

where $G\left(\theta \right)$ is the mean value of leaf projection function for direction $\theta$. It is generally fixed to 0.5 for global-scale remote sensing applications [33]. LAI and CI are the leaf area index and clumping index, respectively. They are directly available from satellite remote sensing products. To explore the uncertainties of the remote sensing products on the estimation of SIF escape probability, artificial errors (−50% and +50% constant value error of LAI*CI) were generated and then were used to estimate SIF escape probability.

2.4. Evalution Process

In Section 3.2, all four estimation formulas in Table 1 were evaluated with ‘real’ values simulated by the DART model (Section 3.1); in Section 3.3, only $\frac{NI{R}_v}{i_{°}·{\omega }_{\lambda }}$ was used to evaluate the usability of remotely sensed canopy interceptance (i_o); in Section 3.4, only $\frac{NI{R}_v}{\mathrm{fAPAR}·{\omega }_{\lambda }}$ was used to evaluate the influence of woody elements on the estimation of SIF escape probability. The complete evaluation process is outlined as a flowchart in Figure 1.

3. Results

3.1. Simulation of ‘Real’ f^esc in Discontinuous Forest Canopies

Using the variable input parameters listed in Table 2, TOC SIF and total emitted SIF were simulated respectively using the DART model. After that, SIF escape probabilities (i.e., ‘real’ f^esc) were calculated according to Equation (2) in Section 2.2.1.

Figure 2 shows the simulated TOC SIF and total emitted SIF, as well as corresponding SIF escape probabilities. It appears that varying input parameters for a given sun direction greatly impacts the simulated TOC SIF, total emitted SIF, and corresponding SIF escape probabilities in all simulations. The variability of canopy structure (i.e., FVC) has a larger impact on TOC SIF (Figure 2a) than total emitted SIF (Figure 2b), which then induced a large variability of SIF escape probabilities from 0.2 to 0.4 (Figure 2c). Even though leaf chlorophyll content has a considerable impact on TOC SIF and total emitted SIF, respectively, SIF escape probabilities remain relatively stable if leaf chlorophyll content varies from 15 to 60 μg/cm² (Figure 2f). On the other hand, with the increase of the solar zenith angle from 20° to 70°, TOC SIF and total emitted SIF decrease at most angles, except from 20° to 30° where the total emitted SIF slightly increases in most cases (Figure 2b,e). At the same time, the corresponding SIF escape probabilities remain relatively stable for all simulations between 30° and 70° (Figure 2c,f).

3.2. Estimation of f^esc and SIF_tot Using Simulated Canopy Interceptance and fAPAR_leaf

i_o and fAPAR_leaf were first simulated using methods described in Section 2.2.2 and Section 2.2.3, and then used to estimate SIF escape probabilities using SIT described in Section 2.1. Figure 3 illustrates the variation of i_o and fAPAR_leaf with the increase of FVC, for Cab = 35 μg/cm², background reflectance = Ref-2. Both i_o and fAPAR_leaf increase with solar zenith angle. It appeared that i_o is usually higher than fAPAR for a solar zenith angle larger than 30°.

Figure 4 shows the relative errors of SIF escape probabilities that are estimated using SIT with DART-simulated values (i.e., ‘real’ values). It appears that FVC generally has a much larger impact on the estimation accuracy of SIF escape probabilities than leaf chlorophyll content. For example, in Figure 4a, varying FVC greatly influences the estimated SIF escape probability with relative errors up to about 180%. In comparison, also in Figure 4a, the relative errors caused by the variation of leaf chlorophyll content were no more than 50%. Furthermore, the varying background spectrum has a considerable impact on the estimation of SIF escape probabilities. With a non-reflecting ground, the relative errors between estimation and simulation were no more than 10%. As for the effect of the solar zenith angle, it is clear that the SIF escape probabilities estimated at high solar zenith angles always have smaller relative errors than those estimated at lower solar zenith angles.

We compared the SIF escape probabilities estimated using canopy interceptance as input and those using fAPAR_leaf as input (i.e., comparison between Figure 4a,c, or Figure 4b,d). No apparent differences were observed in most cases. For the estimations using high FVC (i.e., FVC equal to 0.8), however, the estimation accuracy of SIF escape probability using fAPAR_leaf as input was slightly better than using canopy interceptance as input.

We also compared the estimated SIF escape probabilities using NIR_T as input and using NIR_v as input that was mentioned in Table 1 (i.e., comparison between Figure 4a,b, or Figure 4c,d). The accuracy of the SIF escape probability was better when using NIR_v than using NIR_T. This was especially true with high FVC (i.e., FVC equal to 0.8) as illustrated by the comparison of the estimation accuracy between 0% and 50% with NIR_v (Figure 4b) compared to the accuracy between 50% and 180% with NIR_T (Figure 4a).

We further calculated total emitted SIF ($SI{F}_{tot}^{SIT}$) using the estimated SIF escape probabilities ($f_{SIT}^{esc}$) mentioned above. Here only the results of $\frac{NI{R}_v}{i_{°}·{\omega }_{\lambda }}$ were used to estimate $SI{F}_{tot}^{SIT}$ (corresponding to the $f_{SIT}^{esc}$ in Figure 4b). $SI{F}_{tot}^{SIT}$in Figure 5a was calculated using specific SIF escape probabilities ($f_{SIT}^{esc}$) that were estimated at varying solar zenith angles from 20° to 70°. While for Figure 5b,c, $SI{F}_{tot}^{SIT}$ was calculated using constant SIF escape probabilities ($f_{SIT}^{esc}$) estimated at a fixed solar zenith angle (i.e., SZA = 70° and SZA = 20° for (b) and (c), respectively). In other words, SIF escape probabilities ($f_{SIT}^{esc}$) estimated at SZA = 70° or SZA = 20° were used to assess other cases when solar zenith angles change from 20° to 70°. It is clear that $f_{SIT}^{esc}$ estimated at high solar zenith angles produces higher accuracy (R² = 0.95) of $SI{F}_{tot}^{SIT}$ than that estimated at varying solar zenith angles (R² = 0.85) and low solar zenith angles (R² = 0.81).

3.3. Estimation of f^esc and SIF_tot Using Remotely Sensed Canopy Interceptance

Canopy interceptance (i_o) was estimated using the method described in Section 2.3 with leaf area index (LAI) and clumping index (CI), which were then used to estimate SIF escape probability based on SIT. Figure 6 shows the comparison between estimated and simulated results (i.e., ‘real’ f^esc) of SIF escape probability.

Figure 6a shows the estimated results of SIF escape probability using DART-simulated (i.e., ‘real’) LAI and CI as input for the estimation of canopy interceptance. We found that the estimation accuracy of SIF escape probability was markedly improved (maximum relative error less than 30%) when compared to the estimated results using simulated (i.e., ‘real’) canopy interceptance (maximum relative error up to 50%, see Figure 4b), except for the cases with a non-reflecting background where the SIF escape probabilities were generally underestimated.

Additional errors can be introduced into the estimation of SIF escape probability when remotely sensed data of LAI and CI were used as input compared to real LAI and CI values. To explore this impact, artificial errors (−50% and +50% constant value error of LAI*CI) were generated and then were used to estimate SIF escape probability (Figure 6b,c). Introducing a −50% error of LAI*CI, SIF escape probabilities were overestimated in most simulations, especially for the cases with high FVC and background mostly covered by vegetation (Ref-2). We noted that the introduction of a +50% error on LAI*CI has no apparent impact on the estimation of SIF escape probability, with only a slight underestimation per simulation.

To further evaluate the influence of remotely sensed canopy interceptance on the estimation of SIF escape probability, we calculated $SI{F}_{tot}^{SIT}$ (Figure 7) using the estimated SIF escape probabilities (i.e., $SI{F}_{tot}^{SIT}=SI{F}_{obs}/f_{SIT}^{esc}$) mentioned in Figure 6. When compared with the DART-simulated reference values ($SI{F}_{tot}^{DART}$), additional errors of LAI*CI reduce the estimation accuracy of $SI{F}_{tot}^{SIT}$ (and $f_{SIT}^{esc}$) (R² reduced from 0.91 to 0.66 and 0.88). Further, underestimated LAI*CI produces a larger influence on the estimation accuracy (R² = 0.66) of $SI{F}_{tot}^{SIT}$ (and $f_{SIT}^{esc}$) than overestimated LAI*CI (R² = 0.88).

3.4. Estimation of f^esc Using fAPAR_canopy and fAPAR_leaf

SIF escape probabilities were also estimated using two types of fAPAR (i.e., fAPAR_canopy = fAPAR_leaf+woody and fAPAR_leaf) whose simulation methods were described in Section 2.2.3. Figure 8 shows that the values of SIF escape probability estimated using fAPAR_canopy as input were smaller than that using fAPAR_leaf as input in all simulations. The estimation accuracy of SIF escape probability was improved only for the case with high FVC (FVC equal to 0.8), and the estimation accuracy was decreased for other simulation cases.

4. Discussion

4.1. DART-Based Simulation of f^esc in Discontinuous Forest Canopies

Although several 3-D SIF models can simulate canopy-leaving SIF [16,17,18], SIF escape probability can be simulated directly only by 3-D models that simulate the SIF emitted per unit area (i.e., per leaf). This is especially important for sparse 3-D vegetation scenes [12]. Since it simulates the total and SIF 3-D radiative budgets, we used the DART model (Figure S2) to simulate the total SIF emitted by the canopy simultaneously with canopy-leaving SIF radiance in a sparse 3-D vegetation scenes (see methods in Section 2.2 and results in Section 3.1). Furthermore, we also simulated canopy interceptance (i_o) and fAPAR for the same 3-D scenes (Figure 3), which enabled the calculation of escape probabilities using SIT.

Recently, researchers from [34] developed a relatively simple SIF model based on escape and recollision probability (FluorRTER), which is capable to simulate canopy-leaving and total emitted SIF in 3-D discontinuous canopies. Compared to the DART model, however, the scenes simulated by FluorRTER are much simpler. Tree crowns are cylinders, ellipsoids, or conics, and the heterogeneity within tree crowns is not considered, just like its parent model [35]. This missing complexity is important, as woody elements influence forest reflectance and SIF [25,29,36,37]. With the simulation of fAPAR of leaves and woody elements using the DART model, respectively, we evaluated the impact of woody elements on the simulation and estimation of SIF escape probabilities (see Section 3.4). Although the fine structure of the simulation scene is well described in DART, the computational cost of simulating the 3-D radiative budget of very large regions is too large for practical applications. Then, spatial sampling of 3-D radiative budget should be considered in order to keep reasonable computer time and memory according to the complexity of canopy. Here, we did not meet such limitations because we conducted a theoretical analysis of canopy-leaving and total SIF for assessing the SIF escape probability of a complex but small canopy. It is interesting to note that the very recent bi-directional Monte Carlo mode of DART, called DART-Lux, decreases computer time and memory of simulations of bottom and top of the atmosphere SIF radiance images by a hundredfold [38]. This mode opens new perspectives but must be combined with the standard DART mode because it does not simulate 3-D radiative budgets. Moreover, DART-based simulation can also be used to evaluate and validate other models with simplified canopy structure that are more suitable for global application, such as the FluorRTER model.

4.2. Influence of Background and Solar Zenith Angle on the Estimation of SIF Escape Probability

We compared the simulation values of SIF escape probability (as ‘real’ value) with the estimated values using SIT (Figure 4). The estimation accuracy of SIF escape probability is quite high with a non-reflecting background, both with dense forest (FVC = 0.8) and sparse forest (FVC = 0.3), which means that SIT can be directly used in the estimation of SIF escape probability for forest area with a non-reflecting background. However, in the presence of realistic background instead of the non-reflecting background, the estimation accuracy greatly decreased for each simulation setup, especially for the sparse scene (FVC = 0.3) and for a vegetation-covered background (Ref-3). This is because the reflecting background contaminates the NIR_v used to estimate the SIF escape probability. This result agrees with previous studies [12,22] and the theoretical limit (i.e., black soil problem) of SIT [19,20,39].

Recently, NIR_v (NIR_v = NIR_T × NDVI) was proposed to replace NIR_T for the estimation of SIF escape probability using SIT [12,14], where the use of NDVI is intended to reduce the contamination of soil reflectance on canopy reflectance. Our results confirmed the feasibility of NIR_v for specific simulation cases where the background was not totally covered by vegetation. However, the estimation accuracy of SIF escape probability was not improved with the use of NDVI for the simulation case where the background was mostly covered by vegetation (i.e., Ref-3). Considering the seasonal variation of understory LAI in sparse forests [31,40], the use of NDVI for the estimation of SIF escape probability is still questionable in summer when the background is mostly covered by vegetation (i.e., high understory LAI value).

With regards to the effect of the solar zenith angle on the estimation of SIF escape probability, Figure 4 shows that the accuracies of SIF escape probability estimated at a high solar zenith angle are higher than those estimated at a low solar zenith angle in most studied cases. This phenomenon may be due to the potential link between solar zenith angle and background reflectance. With the solar zenith angle increasing from 20° to 70°, more incident photons will be absorbed by the tree crown, and fewer incident photons will be absorbed by the background, which then will weaken the impact of background reflectance on the estimation of SIF escape probability. The above phenomenon was further confirmed by Figure 5b where the estimation accuracy of total emitted SIF was clearly improved when using a constant value of SIF escape probability estimated at 70° (R² increased from 0.85 to 0.95). Conversely, the estimation accuracy of total emitted SIF slightly decreased (R² decreased from 0.85 to 0.81) when using a constant value of SIF escape probability estimated at 20° (Figure 5c). It is worth noting that these findings were based on the assumption that other input parameters remain constant when the solar zenith angle changes from 20° to 70° (i.e., over the season). In reality, it is not the case, because canopy structure and leaf optical properties change over the season [1,31,41]. Fortunately, Figure 2 shows that the impact of leaf optical properties (represented by Cab) on SIF escape probability is slight in the far-red region, which agrees with results obtained with the SCOPE model (Figure 4 in [22]).

4.3. Influence of Remotely Sensed i_o and fAPAR on the Estimation SIF Escape Probability

In general, remotely sensed canopy interceptance (Figure 6a) led to a relatively accurate SIF escape probability compared to the simulated canopy interceptance (i.e., ‘real’ value, see Figure 4b). Interestingly, the use of remotely sensed canopy interceptance slightly improved the estimation accuracies for configurations with the natural background spectrum. This phenomenon was explained by the joint influences of canopy interceptance and canopy reflectance on the estimation of SIF escape probability, which we term the contamination compensation effect. Firstly, SIF escape probability was generally overestimated by the use of canopy reflectance (NIR_v or NIR_T), because canopy reflectance is contaminated by the reflectance of the background (see Figure 4). Then, the use of overestimated canopy interceptance (Figure S3a–c) partly reduces the estimation errors caused by the contaminated canopy reflectance, especially for the sparse forest canopies. On the contrary, for these configurations with a non-reflecting background spectrum, because there is no more contamination caused by the background reflectance, the estimation accuracy of SIF escape probability was mostly reduced (Figure 6a). However, it should be noted that the above analysis was based on the assumption that $G\left(\theta \right)$ was fixed to 0.5, as for global-scale remote sensing applications [13,14,25,33], and LAI and CI were set to simulation values (i.e., ‘real’ value).

The introduction of artificial constant value errors on LAI*CI has a remarkable influence on the estimation accuracy of SIF escape probability (and total emitted SIF), with a clear overestimate of SIF escape probabilities if LAI*CI is underestimated (Figure 6b). To better explain this trend, canopy interceptance was plotted using LAI*CI with corresponding artificial errors (Figure S3a–c). It demonstrated that an underestimated LAI*CI implies an underestimated canopy interceptance, which amplifies the overestimation induced by the contaminated canopy reflectance. An overestimated LAI*CI tends to produce an overestimated canopy interceptance, which is somewhat compensated by the contaminated canopy reflectance.

The errors of satellite-retrieved LAI and CI products, including theoretical and physical errors, are inevitable [42]. These will be more significant for forest biomes where the accurate estimation of LAI remains challenging with an uncertainty relatively higher than for other biomes [43,44,45]. In a recent study, Zhang et al. [14] demonstrated that total emitted SIF, estimated using satellite-retrieved LAI and CI products, was better correlated with GPP for all C3 sites combined, but was less correlated for forest sites. The above phenomenon could be partly attributed to the high error and uncertainty of satellite-retrieved LAI and CI products. Our result (see Figure 7) also demonstrated the high uncertainty of estimation of canopy interceptance in the forest using LAI and CI, especially for underestimated LAI*CI. Furthermore, $G\left(\theta \right)$ was fixed to 0.5 for all biomes, which can also be a potential source of error because it only holds for homogeneous canopy with spherical leaf orientation distribution [25,33]. Overall, the estimation of canopy interceptance, and the SIF escape probability or total emitted SIF, using satellite-retrieved LAI and CI products is still questionable and needs further investigation.

Owing to the high uncertainty of estimated canopy interceptance using satellite-retrieved LAI and CI products, Zeng et al. [12] proposed to use fAPAR to approximate canopy interceptance instead of canopy interceptance. Our simulation results support this approach since it leads to comparable and even better estimation of SIF escape probability than the method using canopy interceptance (Figure 3). However, fAPAR for whole canopy (i.e., fAPAR_canopy) generally can be divided into fAPAR for photosynthetic elements (i.e., fAPAR_leaf) and fAPAR for non-photosynthetic elements (i.e., fAPAR_woody), respectively [46,47]. The impact of fAPAR_woody is neglected in most satellite-retrieved fAPAR products [48]. As a result, the direct use of satellite-retrieved fAPAR products (mostly fAPAR_canopy) to estimate SIF escape probability still needs further investigation [12,49]. Our simulation results illustrated that using fAPAR_canopy will induce an underestimation of SIF escape probability when compared to fAPAR_leaf (Figure 8). This underestimation can be a potential advantage for the sparse forest where the SIF escape probabilities are generally overestimated due to the contamination of background, providing an additional mechanism behind the contamination compensation effect. Still, more measurement data from stand to global scale are needed to explore this finding and develop a robust method to separate the various canopy layers.

5. Conclusions

This study evaluated the feasibility of SIF escape probability estimation methods in discontinuous forest canopies. In general, these estimation methods based on spectral invariants theory have good accuracy for relatively dense forest canopies. The use of NDVI can partly improve the low estimation accuracy for sparse forest canopies. However, if the forest background is partly or totally covered by vegetation, NDVI failed to eliminate the effects of background reflectance on the estimation of SIF escape probability, because NDVI failed to totally reduce the contamination of background reflectance on canopy reflectance. SIF escape probability estimated at a high solar zenith angle will have a higher accuracy than those estimated at a low solar zenith angle in most studied cases due to the potential link between solar zenith angle and background reflectance. For global-scale remote sensing applications, spectral invariant canopy interceptance can be remotely sensed using the current satellite LAI and CI products. However, additional errors will be introduced with the use of satellite products, especially when the product of LAI and CI was underestimated. Our simulation results also show that canopy interceptance can be replaced by fAPAR with a comparable estimation accuracy, and the direct use of fAPAR_canopy may have a potential advantage than fAPAR_leaf in sparse forest canopies due to a contamination compensation mechanism. To summarize, our results improved our understanding of SIF escape probability estimation methods in discontinuous forest canopies. However, future work is required to evaluate these findings using measurement data. Such work could explore OCO-2 or TROPOMI SIF retrievals or field observations.