Variation in photosynthetic and nonphotosynthetic vegetation along edaphic and compositional gradients in northwestern Amazonia

Introduction Conclusions References


Introduction
Tropical forests vary greatly in aboveground properties such as canopy height, canopy structure, and plant species composition, due to underlying variations in soils and geological formations.Changes in these aboveground properties have been observed over thousands of kilometers (Malhi et al., 2004;ter Steege et al., 2006) and distances as short as 2-3 km (Higgins et al., 2011;Salovaara et al., 2005;Tuomisto et al., 2003b).Due to the inaccessibility of these forests and challenges of sampling in these environments, however, field data on forest properties are often difficult to obtain.Satellite imagery offers a unique opportunity to explore tropical forests for broad-scale variations in aboveground properties that might otherwise be impossible to detect, and to compare these to underlying geological or edaphic properties.
Landsat is the most widely available and commonly used source of satellite imagery for tropical forests, and reflectance in individual Landsat bands is known to correspond to variations in soils and plant species composition (Higgins et al., 2012;Salovaara et al., 2005;Tuomisto et al., 2003a).Analysis of Landsat data can also be used to give insight into the biophysical properties of tropical forests.Individual pixels in Landsat imagery typically represent reflected radiation from a mix of features including live or dead vegetation, bare ground, or human development.Using spectral mixture analysis (SMA), mixed spectra for individual pixels can be decomposed into their component spectra, also known as endmembers (Somers et al., 2011).Typical endmembers for forested ecosystems are photosynthetic vegetation (PV), nonphotosynthetic vegetation (NPV), and bare substrate (Asner et al., 2005b).The products of this process are maps of percentages of PV, NPV, or bare substrate across full Landsat images.
Photosynthetic vegetation endmembers are typically defined from aboveground measurements of intact canopy, collected aerially or from high spatial and spectral resolution orbital sensors such as Earth Observing-1 Hyperion (Asner et al., 2005a).High PV values thus translate into continuous green, leafy canopy with little exposed woody material, litter, or bare ground.Nonphotosynthetic vegetation and bare substrate endmembers are usually calculated using ground-based spectroscopic measurements.Nonphotosynthetic endmembers are defined from ground litter, exposed tree trunk and branches, and deforestation residue; bare substrate endmembers are constructed from a range of soil types with varying water content and organic matter.High NPV thus translates into forests with open canopies and greater amounts of exposed woody or senescent material; high bare substrate values translate into forest removal, exposed soils and rock, or human habitation.
Spectral unmixing has provided valuable insight into the condition of tropical forest vegetation and changes in tropical forests through time.Fractional cover maps of PV, NPV, and bare substrate have been used to identify patterns of selective logging across tropical forests (Asner et al., 2005a;Oliveira et al., 2007;Souza et al., 2005;Allnutt et al., 2013), and to detect subcanopy burn scars (Alencar et al., 2011).In temperate ecosystems, PV values have been found to be superior to the normalized differential vegetation index (NDVI) in remotely sensed calculations of vegetation productivity (Huang et al., 2012).Beyond the tropics, SMA has been used extensively on both global and regional scales for decomposing images into a variety of constituent surface covers (Somers et al., 2011).
Despite the use of fractional cover to detect human disturbance and fire in tropical forests, little is known about natural variations of PV and NPV in intact forests.Recent studies have found broad-scale patterns in plant species composition, soil properties, and underlying geology in western Amazonian forests (Higgins et al., 2011;Tuomisto et al., 2003bTuomisto et al., , 1995;;Salovaara et al., 2005).This raises the possibility that these geological, edaphic, and compositional drivers may translate into changes in the abundance of photosynthetic and nonphotosynthetic vegetation displayed by the canopy, and that this might be detected by spectral unmixing.If true, this would suggest that there are structural and functional differences between these forest types.
Here we tested whether variations in PV and NPV in intact forests in northwestern Amazonia are sensitive to changes in soils and plant species composition.We concentrated our study on two regions of northwestern Amazonia that are known to contain boundaries between two widespread geological formations: the Pebas and Nauta formations.These formations differ substantially in plant species composition and soil properties (Higgins et al., 2011(Higgins et al., , 2012)), making them ideally suited for this purpose.Moreover, we have previously established over 100 plant inventories and soil samples for these study areas, distributed among the two formations.
We used the Carnegie Landsat Analysis System (CLAS) to quantify the fractional cover of PV and NPV throughout both study regions (Asner et al., 2005b), and calculated mean PV and NPV values for all inventories.We then calculated the relationship between PV and NPV values for these inventories, and either soil fertility or plant species composition.We also partitioned the variation in PV and NPV between soils and species composition, in order to determine which was more important in regulating changes in fractional PV and NPV cover.
By comparing these field data with PV and NPV values calculated from Landsat imagery, we calculated the degree to which quantities of photosynthetic and nonphotosynthetic vegetation are controlled by variations in plant species composition, soil properties, or a combination of both.To our knowledge, this is the one of the first times that PV and NPV data have been related to patterns in plant species composition or soil properties in intact tropical forests (see also Carvalho et al., 2013).

Study areas
We conducted our fieldwork in two study areas, Pastaza-Tigre and Curaray, each of which were centered on boundaries between the Pebas and Nauta formations (Fig. 1a-c).The Pebas Formation consists of relatively cation-rich sediments deposited during low-energy semi-marine conditions of the mid-to late Miocene and recently exposed by recent river incision (Hoorn et al., 2010;Räsänen et al., 1995).The Nauta Formation consists of relatively cation-poor sediments deposited by high-energy fluvial conditions during the Plio-Pliocene (Rebata et al., 2006).These two geological formations are widespread across northwestern Amazonia and account for the majority of land surface in northeastern Amazonian Peru (INGEMMET, 2000).
The Pastaza-Tigre study region was located between the Pastaza and Tigre rivers (Fig. 1a and b), and was accessed by a service road for the Lot 1AB oil pipeline.The Curaray study area was located to the southwest of the Curaray River (Fig. 1a and c), and was accessed by helicopter during seismic oil exploration.We established a total of 117 linear transects at these study areas, consisting of 65 transects at Pastaza-Tigre and 52 at Curaray, along which we sampled plants and soils.We restricted our sampling to non-inundated (i.e., terra firme) primary rainforest.Due to the proximity of both study areas to the Equator (between two and three degrees south), the climate at both areas is relatively aseasonal, and these forests do not exhibit pronounced phonology.

PV and NPV data
We used the USGS GLOVIS service (http://glovis.usgs.gov/) to identify and download available cloud and haze-free Landsat images for our study areas, totaling five images for Pastaza-Tigre and three for Curaray (Table 1).These included images from both the Landsat 5 and 7 missions, and for the latter only images prior to the scan line corrector system failure were acquired.We then used the Carnegie Landsat Analysis System (CLAS) to calculate the proportion of each Landsat pixel occupied by PV or NPV for each of the eight Landsat images, based on all six non-thermal bands (Table 1).CLAS uses an automated Monte Carlo unmixing (AutoMCU) approach, in which endmember fractions are calculated from an iterative drawing from spectral endmember libraries until a stable solution is produced (Asner et al., 2005b).CLAS has previously been used to estimate PV and NPV values across the Legal and Peruvian Amazonia, including our study areas, and has been optimized for use in these forests (Asner et al., 2005a;Oliveira et al., 2007).As such, the uncertainties associated with these estimates are low and range from 1 to 4 %, depending on type of vegetation (Asner et al., 2005b).The products of this process were PV and NPV images of the same size and spatial resolution as the original Landsat images, in which pixel values gave the estimated percentage cover of photosynthetic or nonphotosynthetic vegetation.As such, when we refer to photosynthetic or nonphotosynthetic vegetation quantities in this manuscript, we are specifically referring to the percentage cover values estimated by CLAS.In total, we generated PV and NPV images for five dates for the Pastaza-Tigre study area, and three dates for Curaray (Table 1).
We used 250 m buffer areas for each transect (Tuomisto et al., 2003a) to calculate PV and NPV statistics for each transect, and repeated this for all image dates throughout both study regions.Before calculating values for Pastaza-Tigre transects, we edited the buffer areas to remove areas that overlapped with road features.We also omitted, for single image dates, PV and NPV values for any transects that contained clouds or cloud shadows (Table 1).On average, buffer areas for both Pastaza-Tigre and Curaray contained 550 pixels of PV and NPV data, and all analyses were conducted with ArcGIS v. 10 (ESRI Inc., Redlands, CA, USA).

Plant and soil data
The tree species inventory in Amazonian forests is notoriously difficult due to large numbers of species, tall and difficult-to-identify individuals, and poorly known taxonomy (Higgins and Ruokolainen, 2004).Establishing a tree in-ventory network that could sample the full geological and edaphic ranges of our study areas was thus not feasible, and we focused our plant inventories on a single plant group, the pteridophytes (ferns and lycophytes).Pteridophytes are known to capture a majority of the patterns observed in tree inventories at sites across the neotropics (Ruokolainen et al., 1997(Ruokolainen et al., , 2007;;Tuomisto et al., 1995;Jones et al., 2013), and due to their ease of sampling and identification, this method allows an average of one inventory per day.As a result, this group is a common tool for the study of biodiversity patterns in western Amazonia (Higgins et al., 2011;Tuomisto et al., 1995Tuomisto et al., , 2003c)), and has been used in multiple studies to relate remotely sensed imagery to patterns in plant composition in Amazonian forests (Higgins et al., 2012;Salovaara et al., 2005;Tuomisto et al., 2003a, b).
Our field data consisted of 117 linear transects of 5 × 500 m -65 transects at Pastaza-Tigre and 52 at Curaray -along which we had previously collected presence-absence data for pteridophytes (Higgins et al., 2011).Only individuals with at least one leaf (for ferns) or stem (for lycophytes) longer than 10 cm were recorded, and epiphytic and climbing individuals were recorded only if they had green leaves ≤ 2 m above ground.Permits for collection and export of plant specimens for identification were obtained from the Peruvian National Institute of Natural Resources (INRENA), and we deposited vouchers for all species at herbaria in Peru (AMAZ and USM) and Finland (TUR; Thiers, continuously updated, 2014).
In addition, we collected soil samples at 50 m, 250 m, and 450 m along each transect.Each of these three soil samples consisted of five subsamples of the top 10 cm of mineral soil, collected in an area of 4 m × 4 m (Higgins et al., 2011).These five subsamples were located such that one subsample was located in the center of the transect, and the remaining four subsamples were placed at 2 m forward along the transect, back along the transect, to the left , and to the right.These subsamples were then combined in the field into a single sample, and equal dry weights of the three samples were combined into one sample per transect for analysis.Soil samples were analyzed at MTT Agrifood (Jokioinen, Finland) for pH; loss on ignition (a measure of organic matter content); P concentration (Bray method); and extractable Al, Ca, K, Mg and Na (in 1M ammonium acetate).In addition, percentages of sand, silt and clay were determined at MTT Agrifood (Curaray samples; sieving and pipette methods) and the University of Turku Department of Geology (Turku, Finland; laser diffraction).
We used nonmetric multidimensional scaling (NMDS) to reduce our floristic data to a single variable for comparison to soil properties, PV, and NPV values (Higgins et al., 2012).Like other ordination methods, NMDS is a data-reduction technique that identifies the dominant trends in a multidimensional data set -e.g., the simultaneous change in occurrence of multiple plant species along an environmental gradient -and produces a smaller number of dimensions that are more suitable for analyses.In the case of a particularly strong environmental gradient, one NMDS variable can explain the majority of variation between sites in a data set consisting of multiple plant species.We calculated onedimensional NMDS solutions for our plant inventories using the one complement of the Jaccard index as a distance measure.For this analysis we ran a maximum of 400 iterations from 40 random starting configurations, and applied an instability criterion of 10 5 .All NMDS calculations were performed with PC-ORD v. 4.41.

Regression analyses using CLAS data
We used simple linear regressions to calculate the percent of variation in PV and NPV explained by soils and floristic composition for all images in both study regions (Table 1).For these analyses, floristic composition was represented by a single NMDS axis, and cation concentration was represented by the log-transformed sum of four cations (Ca, Mg, Na, and K; abbreviated as LSC).
In addition, due to differences in illumination between image dates, we observed consistently higher PV and NPV values for some dates versus others.To normalize for these between-date variations, and to calculate a single value for comparisons to soil and plant inventory data, we used the difference between PV and NPV (PV -NPV; Asner and Warner, 2003).For each Landsat image date, we subtracted the mean NPV of each transect from its mean PV.We then averaged the PV-NPV values for each transect across all image dates, resulting in a single PV-NPV value for each transect.We then used linear or multiple regressions to calculate the percent of the variation in these PV-NPV values that was explained by NMDS values, LSC values, or the two together (Table 2).These calculations were performed separately for both study areas, and we excluded transects which contained clouds or cloud shadows (as described above).We then used these results to partition the variation in PV-NPV at both study areas into three components (Peres-Neto et al., 2006): variation explained uniquely by variation in LSC values; variation explained uniquely by variation in NMDS values; or variation explained equally by both.We calculated the percentage of variation in PV-NPV explained uniquely by LSC or NMDS values as "C − A" or "C − B", respectively, where "A" represented the percent variation in PV-NPV explained by NMDS values alone (simple linear regression); "B" represented the variation explained by LSC values alone; and "C" represented the variation explained by the combination of LSC and NMDS values (multiple regression).We calculated the percent of variation that was explained equally by LSC and NMDS values as the percentage of variation explained by the combination of LSC and NMDS minus the percentage explained uniquely by both (i.e., "C − (C − A) − (C − B)" simplifying to "A + B − C").

Soil properties and plant species composition
As previously reported (Higgins et al., 2011), the geological patterns in our study area corresponded to abrupt changes in soil cation concentrations.Cation concentrations on the Pebas Formation were 10 times greater than on the Nauta Formation at the Pastaza-Tigre and Curaray region, and 7 times greater at the Curaray region (Fig. 1b and c; sum of concentrations of Ca, Mg, Na, and K).These changes in soils also corresponded to clear changes in plant species composition (Fig. 1b and c).In total, our inventories included 147 and 127 species at Pastaza-Tigre and Curaray, with a mean of 34 and 30 species per transect.Single-axis NMDS ordinations for these pteridophyte inventories explained 80 and 90 % of the floristic pattern in the full data sets at the Curaray and Pastaza-Tigre, respectively, as measured by comparing the distances between inventories based on NMDS scores to the distances between based on the one complement of the Jaccard index.This indicated both strong compositional patterns in these data, and supported the use of these NMDS axes in comparisons to the soil, PV, and NPV data.Linear models using the log-transformed sum of four cations (Ca, Mg, Na, and K) explained 90 and 71 % of the variation in NMDS ordination scores for sites at Pastaza-Tigre and Curaray, respectively (P < 0.001).

Photosynthetic and nonphotosynthetic vegetation
These patterns in geology, soils, and floristic composition were translated into clear patterns of PV and NPV.Forests growing on the Pebas Formation had higher PV and lower NPV values than forests growing on the Nauta Formation (Fig. 1d-g).These differences in fractional cover were strongly explained by differences in cation concentrations between the two geological formations (Fig. 2c, d, g and  h), such that increased cation concentrations resulted in increased PV and decreased NPV.On average, the logtransformed sum of four cations (Ca, Mg, Na, and K) explained 67 and 60 % variation in PV and NPV values at Pastaza-Tigre, and 42 and 65 % at Curaray (P < 0.001 for all regressions; Table 1).
These differences in PV and NPV were equally well explained by variations in plant species composition (Fig. 2a, b, e and f), such that increased NMDS values -corresponding to the transition from the cation-poor Nauta Formation to the richer Pebas Formation -resulted in increased PV and decreased NPV (P < 0.001 for all regressions; Table 1).On average, NMDS scores explained 68 and 59 % variation in PV and NPV values at Pastaza-Tigre, and 41 and 61 % at Curaray.Furthermore, the effects of soil fertility and floristic composition on PV or NPV, as indicated by the slope of the regressions, were identical despite date-to-date variations due to illumination or seasonal effects (Fig. 2a-h).This indicated a strong and consistent relationship between soil fertility, floristic composition, and endmember abundances.Last, due to uncertainties in the estimates of PV and NPV caused by these date-to-date variations in imagery or the CLAS algorithm, these correlations are probably conservative estimates of the strength of these relationships, and the actual correlation values may be greater.
To disentangle the contributions of soil fertility and species composition to PV and NPV results, we combined  PV and NPV into a single metric (PV-NPV), and partitioned the variation in this metric between floristic composition and soil cation concentrations.In combination, soil cation concentrations and NMDS values explained a total of 72 and 67 % of the variation in PV-NPV values at Pastaza-Tigre and Curaray, respectively (Fig. 3, Table 2).Of this, only 3 and 7 % was explained uniquely by soils, and only 1 and 4 % by floristic composition.The large majority (68 and 56 %) of variation in PV-NPV values was explained equally by variation in either soils or floristic composition.This finding suggests that changes in endmember abundance due to cation concentrations were inseparable from changes due to floristic composition (Fig. 3, Table 2).This in turn suggests that variations in PV and NPV due to soil fertility may be mediated by changes in plant species composition.

Discussion
Our analysis of Landsat data for northwestern Amazonia revealed widespread and consistent patterns in photosynthetic and nonphotosynthetic vegetation, corresponding closely to patterns in geology, soils, and plant species composition.Using data from 117 field sites, we found that soil properties and plant species composition explained up to 68 % of variation in PV and NPV in Amazonian forests, and that this increased to 72 % if PV and NPV were combined into a single metric.In addition, the contributions of soil and plant species composition to variations in PV and NPV were inseparable, suggesting that the effect of soils on PV and NPV values may be mediated by variations in species composition.These findings suggest that geology, soils, and plant species composition may produce structural and functional patterns in Amazonian forests.They also suggest that PV and NPV, in  2. addition to their previously demonstrated uses in identifying forest disturbance, might also be used to identify biophysical patterns in Amazonian forests (Alencar et al., 2011;Asner et al., 2005a).
Our findings raise the question of the cause of these patterns and their biological significance.At both study areas, and using eight Landsat image dates, we observed increased NPV and reduced PV when moving from the richer soils of the Pebas Formation to the poorer soils of the Nauta Formation.We believe this may be explained by two factors: thinning of the forest canopy when moving from rich to poor soils, and increased allocation to leaf thickness and defensive compounds.Reduced allocation to leaf area would result in lower PV, and in increased NPV due to increased exposure of woody material such as tree branches or boles.In addition, increased allocation to leaf defensive compounds or thickness might also increase nonphotosynthetic chemical expression in foliage, including structural compounds such as lignin and defensive compounds such as tannins, resulting in increased NPV.These compounds could in turn obscure photosynthetic compounds, resulting in reduced PV.
This scenario is consistent with plot-level studies in Amazonia, which show a trade-off between allocation to plant defensive compounds at low soil cation concentrations; and allocation to light acquisition (e.g., stem and leaf growth) at high soil cation concentrations (Fine et al., 2004(Fine et al., , 2006)).At low cation concentrations plants invest more heavily in defending scarce nutrients, while at high concentrations plants invest more heavily in light acquisition and photosynthesis (Asner et al., 2014).Our findings suggest that these tradeoffs occur at broad scales: from assemblages of species with higher allocations to photosynthetic material on rich-soil formations; to assemblages of species with higher allocations to nonphotosynthetic material (e.g., defensive compounds and thicker leaves) on poor-soil geological formations.
There are two alternative explanations, however, for the observed regional-scale variations in canopy PV and NPV.First, it is possible that gap formation rates on these two geological formations differ such that a higher frequency of gaps on the Nauta Formation translates into higher NPV and lower PV values on this formation.Recent high-resolution lidar studies of canopy structure at these sites, however, indicate that gaps are more frequent on the richer soils of the Pebas Formation (M. A. Higgins and G. P. Asner, unpublished data), suggesting that gap dynamics do not affect PV and NPV values in these forests.
Second, it is possible that terrain differences between the two geological formations are responsible for the observed differences in fractional cover, such that increased slope on one formation causes greater inter-crown shadowing and thus higher NPV and lower PV.Slopes are indeed greater on the Nauta Formation in the Pastaza-Tigre study area, as visible in NASA Shuttle Radar Topography Mission (SRTM) data (Higgins et al., 2011).However, this relationship is reversed at the Curaray study area, where slopes are greater on the Pebas Formation and instead correspond to lower NPV and higher PV.These findings suggest that terrain is not a dominant determinant of PV or NPV values in our study areas.
The patterns we observe might also indicate changes in forest productivity between the two formations.Productivity is calculated from remotely sensed data as the quantity of photosynthetically active radiation absorbed by vegetation (APAR) multiplied by the light-use efficiency of that vegetation (ε; Field et al., 1993;Running, 1990).APAR has traditionally been estimated using the normalized differential vegetation index (NDVI), but PV has recently replaced NDVI in some of these models (Huang et al., 2008(Huang et al., , 2012;;Xiao et al., 2004a, b;Yang et al., 2012).PV is thus directly proportional to productivity with the caveat that PV cannot account for variations in productivity due to changes in light-use efficiency.This suggests that the PV patterns observed here may translate directly into patterns in forest productivity, such that productivity is higher on the richer soils of the Pebas Formation and lower on the poorer soils of the Nauta Formation.
Last, our findings suggest that spectral mixture analysis of multispectral imagery for Amazonian forests may facilitate the detection of more than just changes in canopy cover caused by disturbance and fire.Here we observe clear PV and NPV patterns in intact forests, corresponding to soils and species composition, suggesting that spectral mixture analysis might assist in mapping biodiversity in remote forests.

Figure 1 .
Figure 1.Patterns in geology, soils, plant species composition, and PV and NPV at sites in northern Peru.(A) Geological map for northern Peru (INGEMMET, 2000) with extents of maps in (B) to (G) indicated by solid and dashed lines.Red outline in gray inset indicates position of (A) relative to boundaries of Peru (center) and Ecuador (top left).(B) and (C) Nonmetric multidimensional scaling (NMDS) scores and the log-transformed sum of Ca, Mg, Na, and K concentrations (LSC) overlaid upon Landsat data for Pastaza-Tigre and Curaray study areas, respectively.Points indicate inventory sites and are sized by LSC values and colored by NMDS values, such that larger sizes indicate higher LSC values and blue tones indicate higher NMDS values.Landsat bands four, five and seven correspond to red, green, and blue.Yellow line indicates the boundary between the Pebas and Nauta formations; the Nauta Formation is to the west of the boundary in (B) and inside the boundary in (C).The Pastaza River is located in the far west in (B); the Tigre River runs northwest to southeast through the center of (B); and the Curaray River runs from west to east in the north and east of (C).(D) and (E) PV images for the Pastaza-Tigre and Curaray study areas generated from Landsat imagery for 8 September 2006 and 30 August 2000, respectively.Light tones indicate higher values and dark tones indicate lower values.Extent of (E) is indicated by dashed lines in (A).(F) and (G) NPV images for Pastaza-Tigre and Curaray, as per (D) and (E).

Figure 2 .
Figure 2. Relationship between plant species composition, soil cation concentrations, and productivity at sites in northern Peru.(A) to (H) Linear regressions of the log-transformed sum of Ca, Mg, Na, and K concentrations (LSC) or nonmetric multidimensional scaling (NMDS) values against PV or NPV for Pastaza-Tigre sites (A to D) or Curaray sites (E to H).Each regression line represents a single Landsat image date, and lines are color coded by date.Image dates and correlation coefficients are provided inTable 1. P < 0.001 for all regressions.

Figure 3 .
Figure 3. Contributions of soil fertility and plant species composition to PV-NPV values.Pie charts representing the percentage of variation in PV-NPV values explained by NMDS values alone (red), LSC values (yellow), or equally by both (orange), for the Pastaza-Tigre and Curaray study areas; correlation coefficients are provided in Table2.

Table 1 .
Correlation coefficients (r 2 ) for linear regressions of PV or NPV values against NMDS or LSC, by study area and Landsat image date (all values significant at P < 0.001).Format is YY-MM-DD where YY indicates year, MM indicates month, and DD indicates day.bNumber of cloud or shadow-free transects available for analysis, out of a total of 65 for Pastaza-Tigre and 52 for Curaray. a