# Use of Airborne Hyperspectral Imagery to Map Soil Properties in Tilled Agricultural Fields

Authors: W. Dean Hively; Gregory W. McCarty; James B. Reeves; Megan W. Lang; Robert A. Oesterling; Stephen R. Delwiche
Journal: Applied and Environmental Soil Science
 (2011)

---

## Abstract

Soil hyperspectral reflectance imagery was obtained for six tilled (soil) agricultural fields using an airborne imaging spectrometer (400–2450 nm,~10 nm resolution, 2.5 m spatial resolution). Surface soil samples (n=315) were analyzed for carbon content, particle size distribution, and 15 agronomically important elements (Mehlich-III extraction). When partial least squares (PLS) regression of imagery-derived reflectance spectra was used to predict analyte concentrations, 13 of the 19 analytes were predicted with R2>0.50, including carbon (0.65), aluminum (0.76), iron (0.75), and silt content (0.79). Comparison of 15 spectral math preprocessing treatments showed that a simple first derivative worked well for nearly all analytes. The resulting PLS factors were exported as a vector of coefficients and used to calculate predicted maps of soil properties for each field. Image smoothing with a 3×3 low-pass filter prior to spectral data extraction improved prediction accuracy. The resulting raster maps showed variation associated with topographic factors, indicating the effect of soil redistribution and moisture regime on in-field spatial variability. High-resolution maps of soil analyte concentrations can be used to improve precision environmental management of farmlands.


---

## Body

## 1. Introduction

Spatial assessment of soil properties is important for understanding the dynamics of agricultural ecosystems.  Site specific data can provide information that is critical to maintaining healthy soils and adequate nutrient supply for crop production, preventing losses of nutrients and sediments to the environment, and evaluating the transfer of elements such as carbon between land and atmosphere.  Research has demonstrated that soil properties such as carbon content are correlated with field topography, soil texture, electrical conductivity, and soil reflectance [1–4].  A study by Venteris et al. [5] documented accumulation of carbon in low areas of fields following soil translocation from higher areas, with resulting carbon loss and soil degradation in elevated areas, and Thompson et al. [6] used soil-landscape modeling techniques to evaluate topographic distribution of soil texture and carbon content. These geographic approaches accounted for 28% to 68% of variation in measured carbon and demonstrated the complexity of environmental and management practices that affect soil characteristics.  Recent research into soil health and sustainable cropping systems has demonstrated the potential of improved systems management based on knowledge of distributed soil properties [7].  Contemporary farm management relies on moderate resolution soil maps derived from photo and topographic interpretation. Accurate mapping of soil properties is made difficult due to high spatial variability observed within agricultural fields, errors in spatial assessment of soil properties can result from inadequate or biased sampling of the landscape, and the high cost associated with collecting and analyzing soil samples often limits the amount information available to farmers and land managers.  However, advances in remote sensing technology are now providing tools to support geospatial mapping of soil properties, with applications in agricultural and environmental management.Diffuse reflectance spectroscopy offers a rapid and nondestructive means for measurement of soil properties based on the reflectance spectra of illuminated soil [8–10].  A growing body of literature supports the use of spectral reflectance to determine soil properties, mostly using laboratory instrumentation to measure soil reflectance in the visible (400–700 nm) near infrared (700–2500 nm) and mid-infrared (2500–25,000 nm) wavelengths.  Partial least squares (PLS) regression has emerged as a successful chemometric method for extracting predictive information from spectral reflectance datasets [10–12].  The PLS method, characterizes high leverage orthogonal factors within observed spectral variance and matches them to similar factors that describe observed variance within measurements of a corresponding dependant variable. It has been successfully used to predict the results of soil laboratory analysis for carbon content [13, 14], particle size distribution [13, 15, 16] and elemental nutrient content [14, 17, 18], with results sometimes approaching the analytical accuracy of laboratory tests [12, 19].  A review of 44 studies [20] documented R2 associated with prediction of soil carbon that ranged from 0.45 to 0.98, with a median of 0.86.  Prediction accuracy depends on the signal : noise associated with the spectral data, and, like most analytical calibrations, is also highly influenced by the distribution of values in the measured dataset.  The PLS-derived predictive equations, like most analytical calibrations, are most effective when the unknowns fall with the range of observations used to create the predictive equation, and the best success is obtained when an adequate number of locally obtained samples are included in the calibration data set [12, 21].Advances in sensor technology have enabled satellite and airborne collection of hyperspectral imagery, allowing the acquisition of spectrally detailed geospatial reflectance data at field and landscape scales.  By combining PLS regression of soil properties with reflectance data derived from airborne imagery, high-resolution maps of soil properties can be developed, thus overcoming the inaccuracies associated with geospatial interpolation of soil test data. Reports in the literature, for example, [16, 19, 22], indicate great potential for remote sensing approaches to map surface soil properties.  However, additional research is needed to optimize data analysis procedures and improve prediction ability [19].  Separation of signal from noise is an important part of spectral data processing. Geometric and atmospheric adjustments are first required to derive a geospatially representative map of soil reflectance spectra. The imagery can then be smoothed spatially by averaging adjacent pixels, for example using a 3×3 low-pass filter.  This can reduce the noise that results from random signal variability within the detector array while increasing the signal associated with number of observations. It is also common to smooth the spectra in various ways, often by averaging adjacent wavebands or by calculation of first and second derivatives.  Numerous math pretreatments have been evaluated for application of PLS to spectral reflectance data obtained from agricultural soils [8, 14, 23, 24].In this paper, we have three objectives: (1) to evaluate 30 combinations of spectral math pretreatments and imagery smoothing techniques to identify most effective methods of preparing remote sensing data for partial least squares (PLS) analysis of soil properties (2) to develop and validate PLS predictions of soil concentrations for 19 laboratory analytes based on data extracted from airborne hyperspectral imagery and (3) to export resulting PLS vectors to geospatial imagery processing software and calculate high-resolution raster maps of predicted soil characteristics. Six recently tilled agricultural fields were intensively sampled to provide the calibration data set.


## 2. Materials and Methods

### 2.1. Field Sampling

On 10 April, 2007, we collected 315 soil samples from six fields (Figure1) located on working grain farms on the Eastern Shore of the Chesapeake Bay (Delmarva Peninsula, near Easton, MD).  Each of the fields (Temple 1S, 7.3 ha; Temple 1N, 7.1 ha; Temple 2, 8.9 ha; Temple 3, 18.1 ha; Mason, 14.6 ha; Schrader 9.8 ha) was chosen to provide uniform, smooth, bare-soil conditions, had been recently tilled (moldboard plow, field cultivator, and disk), and had little to no vegetation or plant residue. Soil conditions were moderately dry at the time of sampling, with six days of warm spring weather since the previous substantial rainfall (25 mm on 04 April, 2007). All fields were relatively flat (0% to 5% slope). Soil types included moderately well-drained silt loams (Pineyneck PiA, Mattapex-Butlertown MtA), poorly drained silt loams (Othello Ot), and well-drained sandy loams (Indleside IgB).  Although the majority of Eastern Shore farms are managed using no-till practices, the tilled fields were otherwise typical of regional cash grain crop management strategies.Field sites with soil sampling locations (black points) and SSURGO soil boundaries. See Materials and Methods for soil type designations.
(a)(b)(c)(d)(e)(f)Sampling locations (315 total) were established at approximately 40 m intervals in transects across each field. Transects were established by using a tractor to pull a chisel plow shank through the soil at 20 cm depth. All sampling occurred on 10 April, 2007, during one long day of fieldwork with a large sampling crew.  The soil was already well mixed from moldboard plowing, and this fresh tillage created an area of lightly disturbed soil behind the chisel plow shank from which surface soil samples (~400 g) were collected at each location.  Because the chisel plow shank did not invert the soil, this sampling method approximated conditions at the bare soil surface that were observed by the imaging spectrophotometer. Number of samples ranged from 30 to 86 per field.  Global positioning system (GPS) points were established for each sampling location using a handheld Trimble Geo-XT unit with submeter accuracy, calculated as the average of >20 sequential coordinate readings.Soil samples were air dried (>48 hr) and ground to pass through a 2 mm sieve.  Sand, silt, and clay content were determined using the hydrometer method of particle size analysis [25].  Standard soil nutrient analysis was performed at the University of Delaware Soil Testing Laboratory (http://ag.udel.edu/other_websites/DSTP/), using Mehlich III analysis [26] for potassium (K), calcium (Ca), magnesium (Mg), manganese (Mn), zinc (Zn), iron (Fe), nitrogen (N), phosphorus (P), copper (Cu), boron (B), sulfur (S), and aluminum (Al).  Percent phosphorus saturation was calculated as a ratio of phosphorus to iron and aluminum content. Organic matter content (OM) was determined by loss on ignition. Subsamples were roller-milled for 12 hours prior to elemental analysis for carbon (C) and nitrogen (N) content by dry combustion using a TruSpec CN analyzer (Leco Corp, St. Joesph Mich, USA).  None of the samples contained significant inorganic carbon.  On the same day as field sampling for soils, airborne spectral imagery was acquired.


### 2.2. Imagery Collection

The airborne hyperspectral imaging spectrometer (HyperSpecTIR) used in this study was a push broom sensor developed by the SpecTIR Corporation (Reno Nevada, USA).  It measured irradiance in 178 spectral channels between 400 and 2450 nm with approximately 10-nm resolution [27]. Orthorectification was established to within one pixel (<2.5 m) tolerance.  An upward facing radiation sensor measured incoming solar radiation which was used to calibrate imagery to ground reflectance and limit atmospheric effects to the space beneath the aircraft.The plane was flown at 1800 m altitude with a ground speed of 210 km hr-1, between 10 : 00 and 14 : 00 hrs, on 10 April, 2007. At this altitude, the imagery covered a swath 800 m wide (320 adjacent 2.5 m pixels). Flight paths were flown parallel to the principle plane of the sun. A correction for cross-track illumination was tested but did not improve results and was therefore not adopted. A number of sensor errors were identified where faulty detector elements produced erroneous results within particular wavebands, resulting in along-track striping in the field maps of predicted analyte concentrations.  The 15 most obvious of these errors were corrected with spatial smoothing by replacing each faulty reflectance value with the average value of the two neighboring detector elements within the particular faulty waveband. While it was clear from visual inspection of the predicted imagery that bad detector elements remained, only the 15 most obvious errors were corrected.Soil spectra associated with each sampling location were extracted from the imagery by overlaying GPS point shapefiles (<1 m geopositioning error) of sampling locations and selecting the data associated with each underlying pixel. This was done for both the original imagery (1-pixel extraction) and spatially smoothed imagery to which a3×3 low-pass filter had been applied (9-pixel extraction).  Imagery processing was conducted using ENVI 4.7. Elevation data for each pixel were derived from a 3-m LIDAR digital elevation map that was resampled at 2.5 m resolution, with a vertical accuracy of 0.20 m.  A 2.5 m resolution wetness index was also calculated from the LIDAR data, using SAGA software to apply two consecutive iterations of an enhanced lee 3×3 filter.  This measure approximates the cumulative influence of upslope contributing area. Near infrared reflectance spectra for each of the 315 dried, ground soil samples were also obtained in the laboratory using a bench spectrometer and controlled light source, with methods and PLS results reported in McCarty et al. [28].


### 2.3. Spectral Data Processing

Fifteen spectral math pretreatments were evaluated using a SAS Ver. 9.12 program modified for hyperspectral data processing [23].  These included untransformed spectra, first derivatives with gap ranging from 1 to 64, and second derivatives with gap ranging from 1 to 64. Each of these math pretreatments was applied separately to the smoothed and unsmoothed spectra, resulting in a total of 30 data combinations that were used independently to calculate PLS predictions for each laboratory analyte using all 315 samples. Both the spectral data and the analyte values were mean centered prior to PLS analysis. The number of factors used in each PLS regression was chosen by the PRESS algorithm [24] within SAS Proc PLS and ranged from four to ten.  Testing for outliers was not performed, and all observed values were included in the analysis. Goodness of fit was determined using repeated leave-one-out cross validation, with results presented as coefficient of determination (R2). Ranking analysis was employed to test for significant differences among math treatments, for all analytes predicted with R2>0.50.Once the best math treatment was chosen, PLS Toolbox Ver. 4.0 (Eigenvector Research, Wenatchee, Wash, USA), operating within the Matlab (Ver. 7.0) environment, was used to perform PLS analysis using mean-centered spectral and analyte data from five of the six fields (269 samples). Goodness of fit was determined using repeated leave-one-out cross validation, with results presented as coefficient of determination (R2) and residual mean square error (RMSE).  Data from the remaining field (Temple 1S, 46 samples) were treated as unknown samples and predicted from PLS coefficients, with prediction accuracy reported as bias (mean of predicted values minus mean of observed values) divided by the mean of predicted values (bias/mean(pred)) as well as standard deviation of validation [24, Section 18.8] divided by standard error of prediction (sd/se pred).The use of Matlab allowed the mathematical flattening of the PLS factors into a 178 band vector of coefficients, and the scalar product of this vector with each pixel’s imagery-derived spectrum was used to calculate predicted analyte concentrations. This calculation was made for each pixel of the hyperspectral imagery by using IDL code run within the ENVI 4.7 programming environment to calculate the appropriate band math. In this manner, geospatial field maps of predicted analyte concentrations were produced.


## 2.1. Field Sampling

On 10 April, 2007, we collected 315 soil samples from six fields (Figure1) located on working grain farms on the Eastern Shore of the Chesapeake Bay (Delmarva Peninsula, near Easton, MD).  Each of the fields (Temple 1S, 7.3 ha; Temple 1N, 7.1 ha; Temple 2, 8.9 ha; Temple 3, 18.1 ha; Mason, 14.6 ha; Schrader 9.8 ha) was chosen to provide uniform, smooth, bare-soil conditions, had been recently tilled (moldboard plow, field cultivator, and disk), and had little to no vegetation or plant residue. Soil conditions were moderately dry at the time of sampling, with six days of warm spring weather since the previous substantial rainfall (25 mm on 04 April, 2007). All fields were relatively flat (0% to 5% slope). Soil types included moderately well-drained silt loams (Pineyneck PiA, Mattapex-Butlertown MtA), poorly drained silt loams (Othello Ot), and well-drained sandy loams (Indleside IgB).  Although the majority of Eastern Shore farms are managed using no-till practices, the tilled fields were otherwise typical of regional cash grain crop management strategies.Field sites with soil sampling locations (black points) and SSURGO soil boundaries. See Materials and Methods for soil type designations.
(a)(b)(c)(d)(e)(f)Sampling locations (315 total) were established at approximately 40 m intervals in transects across each field. Transects were established by using a tractor to pull a chisel plow shank through the soil at 20 cm depth. All sampling occurred on 10 April, 2007, during one long day of fieldwork with a large sampling crew.  The soil was already well mixed from moldboard plowing, and this fresh tillage created an area of lightly disturbed soil behind the chisel plow shank from which surface soil samples (~400 g) were collected at each location.  Because the chisel plow shank did not invert the soil, this sampling method approximated conditions at the bare soil surface that were observed by the imaging spectrophotometer. Number of samples ranged from 30 to 86 per field.  Global positioning system (GPS) points were established for each sampling location using a handheld Trimble Geo-XT unit with submeter accuracy, calculated as the average of >20 sequential coordinate readings.Soil samples were air dried (>48 hr) and ground to pass through a 2 mm sieve.  Sand, silt, and clay content were determined using the hydrometer method of particle size analysis [25].  Standard soil nutrient analysis was performed at the University of Delaware Soil Testing Laboratory (http://ag.udel.edu/other_websites/DSTP/), using Mehlich III analysis [26] for potassium (K), calcium (Ca), magnesium (Mg), manganese (Mn), zinc (Zn), iron (Fe), nitrogen (N), phosphorus (P), copper (Cu), boron (B), sulfur (S), and aluminum (Al).  Percent phosphorus saturation was calculated as a ratio of phosphorus to iron and aluminum content. Organic matter content (OM) was determined by loss on ignition. Subsamples were roller-milled for 12 hours prior to elemental analysis for carbon (C) and nitrogen (N) content by dry combustion using a TruSpec CN analyzer (Leco Corp, St. Joesph Mich, USA).  None of the samples contained significant inorganic carbon.  On the same day as field sampling for soils, airborne spectral imagery was acquired.


## 2.2. Imagery Collection

The airborne hyperspectral imaging spectrometer (HyperSpecTIR) used in this study was a push broom sensor developed by the SpecTIR Corporation (Reno Nevada, USA).  It measured irradiance in 178 spectral channels between 400 and 2450 nm with approximately 10-nm resolution [27]. Orthorectification was established to within one pixel (<2.5 m) tolerance.  An upward facing radiation sensor measured incoming solar radiation which was used to calibrate imagery to ground reflectance and limit atmospheric effects to the space beneath the aircraft.The plane was flown at 1800 m altitude with a ground speed of 210 km hr-1, between 10 : 00 and 14 : 00 hrs, on 10 April, 2007. At this altitude, the imagery covered a swath 800 m wide (320 adjacent 2.5 m pixels). Flight paths were flown parallel to the principle plane of the sun. A correction for cross-track illumination was tested but did not improve results and was therefore not adopted. A number of sensor errors were identified where faulty detector elements produced erroneous results within particular wavebands, resulting in along-track striping in the field maps of predicted analyte concentrations.  The 15 most obvious of these errors were corrected with spatial smoothing by replacing each faulty reflectance value with the average value of the two neighboring detector elements within the particular faulty waveband. While it was clear from visual inspection of the predicted imagery that bad detector elements remained, only the 15 most obvious errors were corrected.Soil spectra associated with each sampling location were extracted from the imagery by overlaying GPS point shapefiles (<1 m geopositioning error) of sampling locations and selecting the data associated with each underlying pixel. This was done for both the original imagery (1-pixel extraction) and spatially smoothed imagery to which a3×3 low-pass filter had been applied (9-pixel extraction).  Imagery processing was conducted using ENVI 4.7. Elevation data for each pixel were derived from a 3-m LIDAR digital elevation map that was resampled at 2.5 m resolution, with a vertical accuracy of 0.20 m.  A 2.5 m resolution wetness index was also calculated from the LIDAR data, using SAGA software to apply two consecutive iterations of an enhanced lee 3×3 filter.  This measure approximates the cumulative influence of upslope contributing area. Near infrared reflectance spectra for each of the 315 dried, ground soil samples were also obtained in the laboratory using a bench spectrometer and controlled light source, with methods and PLS results reported in McCarty et al. [28].


## 2.3. Spectral Data Processing

Fifteen spectral math pretreatments were evaluated using a SAS Ver. 9.12 program modified for hyperspectral data processing [23].  These included untransformed spectra, first derivatives with gap ranging from 1 to 64, and second derivatives with gap ranging from 1 to 64. Each of these math pretreatments was applied separately to the smoothed and unsmoothed spectra, resulting in a total of 30 data combinations that were used independently to calculate PLS predictions for each laboratory analyte using all 315 samples. Both the spectral data and the analyte values were mean centered prior to PLS analysis. The number of factors used in each PLS regression was chosen by the PRESS algorithm [24] within SAS Proc PLS and ranged from four to ten.  Testing for outliers was not performed, and all observed values were included in the analysis. Goodness of fit was determined using repeated leave-one-out cross validation, with results presented as coefficient of determination (R2). Ranking analysis was employed to test for significant differences among math treatments, for all analytes predicted with R2>0.50.Once the best math treatment was chosen, PLS Toolbox Ver. 4.0 (Eigenvector Research, Wenatchee, Wash, USA), operating within the Matlab (Ver. 7.0) environment, was used to perform PLS analysis using mean-centered spectral and analyte data from five of the six fields (269 samples). Goodness of fit was determined using repeated leave-one-out cross validation, with results presented as coefficient of determination (R2) and residual mean square error (RMSE).  Data from the remaining field (Temple 1S, 46 samples) were treated as unknown samples and predicted from PLS coefficients, with prediction accuracy reported as bias (mean of predicted values minus mean of observed values) divided by the mean of predicted values (bias/mean(pred)) as well as standard deviation of validation [24, Section 18.8] divided by standard error of prediction (sd/se pred).The use of Matlab allowed the mathematical flattening of the PLS factors into a 178 band vector of coefficients, and the scalar product of this vector with each pixel’s imagery-derived spectrum was used to calculate predicted analyte concentrations. This calculation was made for each pixel of the hyperspectral imagery by using IDL code run within the ENVI 4.7 programming environment to calculate the appropriate band math. In this manner, geospatial field maps of predicted analyte concentrations were produced.


## 3. Results and Discussion

### 3.1. Soil Test Results

Observed analyte concentrations for the six sampled fields are described in Table1. Soils were dry at the time of sampling, with moisture content ranging from 10% to 23% for the majority (91%) of samples (Table 1).  Overall, the distribution of observed carbon content values was not large, ranging from 0.6 to 2.0% (Table 1).  The Mason site, which had recently transitioned from long term conventional grain production to organic grain production, exhibited low soil carbon (mean 0.7%) content relative to the other fields (means of 1.2 to 1.4%).  The Schrader site had received long-term applications of dairy manure and exhibited somewhat increased soil C (mean 1.4%) relative to the other fields. The remaining four fields (Table 1: Temple 1S, Temple 1N, Temple 2, and Temple 3) were managed by one farmer, and were under similar management (conventional corn-wheat/soybean rotation). Overall, carbon content was low, as is typical of Maryland Eastern Shore farmland. Concentrations of the remaining analytes were within the normal range for agricultural soils, although variability among sites was not great.Table 1
Observed laboratory analyte concentrations for the six bare soil fields that were sampled on 10 April, 2007, with number of samples collected(n). See Section 2 for analyte descriptions.
FieldC%Sand%Silt%Clay%pHOM%K mg/kgCa mg/kgMg mg/kgMn mg/kgZn mg/kgFe mg/kgAl mg/kgN mg/kgP mg/kgCu mg/kgB mg/kgS mg/kgAl mg/kgMoisture%Elevation mWetnessTemplemin1.0736.528.93.34.61.70943785581.1973270.03231.80.2120.53485−5.07.21Smax1.6365.451.812.56.42.7097412142421116.94097980.125856.00.7982.0798352.712.7(n=46)Mean1.3248.642.49.05.72.22247817142342.52055010.07953.30.4041.4503150.010.27.3 hasd0.157.36.32.00.40.2515417139311.1751210.021100.80.1412.511941.51.1Templemin0.8932.324.65.84.71.30835016991.31152990.02221.90.1927.729911−4.87.11Nmax1.7566.955.716.66.33.0049312522571153.93609600.121755.00.7475.0960202.713.3(n=52)mean1.2149.540.99.65.32.12226785137322.52075590.06673.30.3746.3557160.010.47.1 hasd0.177.66.72.20.40.338919246290.7711460.02290.60.1413.914821.51.2Templemin0.8834.124.73.95.41.404348210740.81293060.03111.90.1726.630611−5.06.82max1.9168.255.417.37.03.105711423420224.44036480.112655.90.67108.8648421.512.7(n=55)mean1.2947.642.59.86.32.17164953256101.92314630.07493.30.4148.2463170.09.88.9 hasd0.227.66.22.40.40.351031856330.774680.02390.90.1319.27051.21.0Templemin1.0019.031.46.35.61.507161513481.4913520.03222.00.3523.83528−5.95.73max2.0459.567.919.77.23.4044518574881128.738010490.112439.14.42147.31049341.012.4(n=86)mean1.3232.954.912.26.32.50210992226343.01836520.07954.00.9337.1652140.010.418.1 hasd0.188.07.12.90.30.377027164231.2771880.01511.10.7714.518840.51.0Masonmin0.5547.69.51.35.70.504335775202.5715670.02342.20.2917.15676−2.76.4(n=30)max1.4087.743.410.96.82.202137752239412.213110330.0732021.20.6129.81033182.913.414.6 hamean0.7569.423.76.96.30.97111540133515.71028030.031484.40.4221.7803100.010.5sd0.199.77.92.50.30.394211134201.9171250.01713.40.083.412531.01.3Schradermin1.0937.218.06.04.92.1025347899224.11205410.04863.50.4428.15416−2.77.9(n=46)max1.8873.945.918.76.83.5060010662121118.64108220.0923212.90.8544.1822232.911.99.8 hamean1.4151.937.710.46.12.65449801160526.01776690.071336.60.6235.6669150.010.6sd0.188.66.72.70.40.348911327191.062660.01351.80.103.46631.00.5Overallmin0.5519.09.51.34.60.50433575540.8712990.02111.80.1717.12995−5.95.7(n=315)max2.0487.767.919.77.23.50974185748811512.241010490.1258521.24.42147.31049422.913.465.8 hamean1.2646.843.110.16.02.21236855186343.41905980.06934.10.5739.5598150.010.3sd0.2513.111.22.90.50.5613623870261.8771680.01531.30.3012.312741.11.0


### 3.2. Choice of Math Treatment

A SAS “shotgun approach” for hyperspectral data processing has previously been used to provide a factorial comparison of the effects of various math pretreatments on PLS analysis of soil spectra [23].  This method was used to compare a total of 15 different spectral pretreatments including 1st and 2nd derivatives with various gap widths (Table 2). Results showed that there were no statistical differences among the majority of the treatments, with the simplest treatments (no derivative, first derivative gap 2) often resulting in the best fit.  This led to the conclusion that the PLS data mining techniques are capable of extracting the majority of signal information from the untransformed spectral data, without the need for math pretreatment.Table 2
Partial least squares (PLS) prediction model goodness of fit (R2) associated with each of 15 math treatments, for the 13 analytes that predicted with R2>0.5, calculated using data from all 315 soil sampling locations1. The first derivative gap two (1Dg2, depicted in bold) was selected as the overall most successful model although it was occasionally outperformed (italic). See Section 2 for analyte descriptions.
DerivativeGapCSandSiltClaypHOMKCaMgMnZnFeAlaverageR2R2R2R2R2R2R2R2R2R2R2R2R2R2NON00.5780.7620.7630.5850.4420.6850.5550.6690.6920.6420.6700.7540.7770.6591ST10.5550.7700.7610.5650.5170.7060.5220.6540.7190.6240.6060.7070.7990.6541ST20.5910.7630.7630.6170.5490.7170.5140.6760.7080.6380.6470.7400.7820.6701ST40.5950.7540.7400.5960.4510.6920.5130.6300.6160.5650.6750.7030.7730.6391ST80.5840.7480.7390.6000.4130.6680.5780.6350.6420.5800.6360.6370.7270.6301ST160.5880.7400.7440.5500.4420.6720.5470.6100.6190.5670.5910.6660.7250.6201ST320.5250.7300.7220.5540.3810.6400.5200.6480.6190.5420.5800.6690.7650.6071ST640.3380.6230.6190.4700.3170.5370.4130.5190.5800.4880.5680.6360.6690.5212ND10.5420.5850.5970.4240.3610.5260.4840.7080.7050.5290.5990.6650.7290.5732ND20.5070.6810.6820.3720.4990.5660.4930.6650.6640.5660.5910.6730.7690.5952ND40.5760.7140.7240.4090.4490.6600.4840.6720.6690.6060.6130.7490.7520.6212ND80.5360.6980.6890.5330.4720.6020.5280.5950.6310.5420.6380.7430.7530.6122ND160.5800.7290.7250.5730.4870.6610.4980.5900.6400.5980.6170.6800.7120.6222ND320.5640.7150.6980.6090.4320.6310.5480.6220.5920.5480.6060.6390.6880.6072ND640.5250.5960.5890.4210.1560.5990.5110.5560.4640.5010.5380.6020.6140.5131Analytes that predicted poorly R2  (<0.5) included: N (<0.303), P (<.355), Cu (<0.358), B (<0.169), S (<0.282), and P saturation (<0.127).There is some argument to be made that different math treatments are more appropriate for particular analytes, due to the physical interaction of light with those particular constituents.  Some evidence for this was shown for potassium, for which a first derivative gap 8 increasedR2 from 0.514 (first derivative gap 2) to 0.578 (Table 2).  However, caution must be taken to prevent overfitting of datasets, and further study is needed to justify the selection of a diversity of math treatments for use with particular analytes.  Ultimately, the first derivative gap 2 was selected as the best math overall treatment even though it was slightly, but not significantly, outperformed in several cases (Table 2).


### 3.3. Partial Least Squares Regression: Calibration Dataset

Once the best math treatment was decided upon, PLS regression was performed in Matlab using a calibration dataset of 269 samples (all sampling locations from five of the six fields), with results reported in Table3. It should be noted that spectral data for the 315 samples used to calculate the SAS PLS in Table 2 were extracted from an earlier version of image processing output in which the 15 band errors had not yet been corrected and to which cross-track illumination had been applied. This discrepancy resulted in slightly lower R2 being used in Table 2 than in Table 3, but was not thought to have affected the relative performance of math treatments. Using Matlab PLS on 269 samples, 13 of the 19 analytes were predicted with R2>0.50 (Table 3(a)), and the remaining seven were predicted with R2<0.40: N (<0.30), P (<0.36), Cu (<0.36), B (<0.17), S (<0.28), and P saturation (<0.13).  A comparison of observed and predicted values for a selection of analytes (carbon, silt, aluminum, and iron) is shown in Figure 2.Table 3
Partial least squares (PLS) model accuracy in predicting soil analyte concentrations1 for (a) the 269 calibration samples using repeated leave-one-out cross validation and (b) 46 validation samples from the field (Temple 1S) that had been left apart from the calibration. Results were derived using a first derivative gap 2 math pretreatment, using spectra derived from unsmoothed imagery (2.5 m2 pixel size) or from spatially smoothed imagery to which a 3×3 low pass filter had been applied. See Section 2 for analyte descriptions. Units refer to analyte residual mean squared error (RMSE) values. Bold indicates sd/sepred >1.20.
Factors no.C%Sand%Silt%Clay%pHOM%K mg/kgCa mg/kgMg mg/kgMn mg/kgZn mg/kgFe mg/kgAl mg/kg8910910104891081010(a) 269 calibration samples collected from five agricultural fieldsusing unsmoothed imageryR20.650.790.790.660.510.750.590.690.690.620.640.750.76RMSE0.197.96.92.20.40.489.5166.150.319.61.449.3104.7using smoothed imageryR20.640.800.800.680.580.770.610.710.710.670.670.780.81RMSE0.187.26.22.00.40.386.1151.845.517.81.243.889.0(b) 46 validation samples collected from the remaining agricultural field (Temple 1S)using unsmoothed imagery2se pred0.209.07.32.40.50.4147.1178.151.427.51.252.9123.93bias0.10−2.881.421.000.030.3424.9240.4531.843.230.452.12107.33bias/mean (pred)0.08−0.060.030.110.010.150.100.050.220.090.180.010.214sd/se pred0.740.810.860.840.890.631.050.960.761.150.921.410.98using smoothed imageryse pred0.159.57.92.90.40.3154.9233.658.538.31.398.0170.4bias−0.10−0.87−0.660.95−0.25−0.13−23.06−2.46−2.07−0.380.10−5.2677.66bias/mean (pred)−0.08−0.02−0.020.10−0.05−0.06−0.100.00−0.02−0.010.04−0.030.14sd/se pred0.861.291.261.070.970.741.041.341.411.250.981.911.291Analytes that predicted poorly (R2 < 0.5) in the set of 315 samples (N, P, Cu,  B, S, and P saturation) are not included here2se pred: sd(observed-predicted)3bias: mean(pred)-mean(obs)4sd: standard deviation.Comparison of observed analyte concentrations with values predicted using partial least squares (PLS) regression on reflectance spectra extracted from smoothed imagery to which a3×3 low-pass filter had been applied, for (a) soil carbon, (b) silt, (c) iron, and (d) aluminum content.  The field containing 41 sampling points that were left out for validation (Temple 1S) is depicted with solid circles.
(a)(b)(c)(d)Prediction accuracy for carbon (R2=0.65) fell well within the range of results (0.45 to 0.98) found in a survey of 44 studies [20] but was somewhat poor in comparison with results found in some other studies for example, [11, 14, 15].  The somewhat poor R2 for carbon might lead to the conclusion that the sensor did not capture a good reflectance signal due to the increased noise often associated with airborne sensors that derives from atmospheric effects and variations in sensor-soil-sun geometry across the imagery. Indeed, a number of studies have found remote sensing spectroscopy to have reduced signal : noise relative to laboratory-based measurements [24].  However, the carbon PLS results reported in McCarty et al. [28], calculated for all 315 soil samples in this data set, using repeated leave-one-out correlation to determine goodness of fit, showed that spectra from the airborne sensor were as effective (R2=0.67) as spectra from a near-infrared benchtop laboratory spectrophotometer (R2=0.64) in predicting soil carbon concentrations.  For the remaining 12 analytes under consideration, the airborne sensor exhibited decreased accuracy relative to the laboratory spectrophotometer in three cases (change in R2 of −0.03 to −0.08) and exhibited increased accuracy in nine cases (change in R2 of 0.01 to 0.20).These observations led to the conclusion that the airborne sensor provides a viable option for mapping soil properties and that the somewhat poor prediction accuracies observed in this experiment apparently stemmed from features associated with the local soil environment or the calibration dataset, rather than the effectiveness of the airborne sensor. The distribution and range of observed analyte concentrations within a calibration dataset can have a substantial impact on prediction accuracy, and the somewhat low prediction accuracy for carbon content observed in this study (R2=0.65) is likely associated with the limited range and low values of soil carbon contents found within the calibration dataset (0.6% to 2.0%).For the other analytes, prediction accuracies (Table3(a)) were comparable with those found in other studies [8, 15, 17, 18]. Using principal components regression, Chang et al. [18] predicted Melich III analyte concentrations, and, similarly to this study, found Cu and P to have poor predictions, silt and clay content to predict with R2~0.8, and a number of other analytes, including Fe, Mg, Mn, K, and pH, predicting with R2>0.6.  It is not always known whether successful predictions are the direct action of the analyte upon the reflectance signal or instead the covariation of analyte concentrations with other factors which are influencing the spectral response [29]. Phosphorus, for example, has no expected reflectance resonance and is generally poorly predicted (e.g., 0.36 in this study, 0.40 in [18], and as low as 0.10 in other studies [8, 29]), and yet occasionally will be predicted with considerable accuracy [8, 29], likely due to covariation with spectrally responsive factors associated with labile organic matter.


### 3.4. Validation

Chemometric predictions of soil properties are typically validated by predicting analyte concentrations for samples that were not included in the calibration data set, and comparing the predicted results to observed values.  This can be achieved by calculating repeated leave-one-out or leave-ten-out predictions, or by selecting a distinct set of samples either randomly, spatially, or based on even sample distribution within the range of observed values [24, Section 18].  We chose to adopt a rigorous validation by removing an entire agricultural field comprising 46 samples (Temple 1S) from the data set, leaving 269 samples for calibration.  This field was one of four Temple fields that received similar crop management (corn-soy/wheat rotation with full tillage and no use of cover crops) and it exhibited analyte concentrations in the mid-range of the six fields (Table 1, Figure 2). When predicted values were calculated for the 46 validation samples and compared to observed values, bias ranged from 1% to 22% of observed analyte concentrations (Table 3(b)).  Seven of the analytes exhibited sd/se pred >1.20 (bold text in Table 3(b)) using the spectra extracted from smoothed imagery, indicating that the prediction was useful for those elements (sand, silt, Ca, Mg, Mn, Fe, and Al), while only one analyte (Fe) exhibited sd/se pred > 1.20  using the spectra extracted from unsmoothed imagery, indicating that spatial smoothing helped to increase signal : noise in the spectral data set.


### 3.5. Prediction Maps

For each analyte, a PLS regression vector of 178 coefficients was exported from Matlab and applied to the mean-centered first derivative reflectance spectra associated with each pixel of the hyperspectral imagery to calculate a 2.5-m raster map of predicted analyte concentrations. Prediction maps for a selection of analytes (carbon, silt, iron, and aluminum) were derived in this manner using both unsmoothed and smoothed imagery (Figure3). Ben-Dor et al. [19] have discussed the need for attention to imagery processing methods to extract maximum information from spectroscopic imagery. In this study, smoothing the hyperspectral imagery with a 3×3 low band filter prior to spectral data extraction resulted in somewhat improved prediction accuracies in comparison to unsmoothed imagery (Table 3(a)).  While this result implies that the smoothed extent of nine adjacent pixels (56.2 m2) provided a better average representation of soil characteristics at each sampling point than did individual 6.25 m2 pixels, it is more likely that the increased accuracy associated with the smoothed imagery is attributable to improvement of signal and reduction in sensor detector element variability (noise) that results from coaddition of the nine adjacent spectra.  Smoothing also generally improved the range and distribution of predicted analyte values found within each field image (Table 4), likely due to the correction of aberrations within the detector elements. However, it should be noted that discontinuities (striping) along the plant’s flight path were evident in both smoothed and unsmoothed predicted imagery (Figure 3), indicating the effects of variable detector element sensitivity along the sensor array and the need for continued improvement in sensor technology.  Several additional imagery adjustments were considered, including cross-track illumination, transformation from reflectance (R) to log(1/R), and also spectral smoothing (not shown). However, none of these operations improved results and they were therefore avoided.Table 4
Predicted concentrations of select analytes (C, Silt, Al, Fe)  derived fromn pixels of near infrared hyperspectral imagery for each of the six bare soil fields using either unsmoothed (1-pixel) imagery with 2.5 m resolution, or spatially smoothed (9-pixel) imagery to which a 3×3 low pass filter had been applied1. See Section 2 for analyte descriptions.
CSiltAlFeField1-pixel%9-pixel%1-pixel%9-pixel%1-pixel mg/kg9-pixel mg/kg1-pixel mg/kg9-pixel mg/kgTemple 1Smin0.971.1314.614.4342454056n=11731max1.861.9480.158.71017947474394mean1.411.4842.141.1633631194192sd0.100.105.95.098924644Temple 1Nmin0.850.961.05.6−311473764n=11436max1.952.0468.959.11304897452429mean1.231.2543.642.1589555200210sd0.110.115.84.91151126065Temple 2min−0.480.8114.68.734−3894054n=14221max2.072.9280.175.710171013474529mean1.271.2842.138.2633512194226sd0.250.185.97.1981464663Temple 3min0.390.600.21.2−1370−288−6n=28955max2.292.4880.874.911731031465474mean1.271.2752.351.7632639173170sd0.170.167.97.51541457163Masonmin0.410.44−7.6−7.5−19829000n=23351max2.271.9464.948.512311180391391mean0.810.8127.728.7755776111108sd0.110.106.47.1961003227Schradermin0.910.9812.313.42552965866n=15684max2.082.0561.656.9996971346315mean1.451.4838.638.3669681177172sd0.130.156.05.887813534Overallmin0.971.1314.614.42552965866n=105378max1.861.9461.648.5996897346315mean1.201.2241.340.56606481681711See Table 3 for the accuracy (R2) associated with each analyte prediction model.Maps of predicted soil carbon content calculated from (a) unsmoothed imagery (1-pixel data extraction) and (b) spatially smoothed imagery (9-pixel data extraction).  Two of the six sampled fields are depicted here, with Temple 3 on the left and Schrader on the right. Predicted values ranged from 0.4% (black) to 2.5% (white).
(a)(b)


### 3.6. Topographic Analysis

Visual interpretation of field topographic features relative to the predicted analyte maps revealed accumulation of carbon, silt, iron, and aluminum in low areas (Figure4).  This result is in agreement with previous studies in Iowa that have linked increased soil C carbon contents in low areas with soil redistribution processes [3], and with a study by Terra et al. [1] who detected a correlation of soil properties with topographic indices and electrical conductivity that explained 50% of observed variability in carbon content.  Other studies have evaluated the links between topography and soil constituents such as organic carbon at the landscape scale, showing accumulation at the bottom of slopes and in valley bottoms [30].  This study’s work was conducted at the field scale, reflecting microtopographic soil distribution processes within a relatively flat elevation gradient (0%–5% slope), as well as the formative effects of soil moisture distribution and spatially variable biological processes resulting from soil heterogeneity and the balance of redoxomorphic status between oxic and anoxic conditions. Understanding distribution patterns of soil analytes, and the influence of hydrogeomorphic controls on these patterns, can provide greater certainty about the influence of soil erosion and soil ecology on the fate of carbon and nutrients.Map of predicted values for select analytes (C, Silt, Fe, Al), overlaid on a high-resolution digital elevation map, for one of the six sampled fields (Mason). Predicted values were derived from imagery that had been smoothed with a3×3 low pass filter.
(a)(b)(c)(d)To investigate potential causes of spatial distribution in predicted analyte values, a 2.5 m elevation map (re-sampled from 3-m LIDAR DEM) was used to calculate a normalized elevation value (observed—field mean elevation) and a soil wetness index (two consecutive iterations of a SAGA enhanced lee3×3 filter) for each pixel of the predicted analyte maps. The normalized elevation variable was somewhat correlated (Table 5(a)) with predicted concentrations of aluminum (R2=0.23, 0.39 for unsmoothed and spatially smoothed imagery, respectively) and was not well correlated with predicted concentrations of silt, iron, and carbon.  The wetness index was similarly correlated with predicted aluminum concentration (R2=0.16,0.32) and was a better predictor for soil carbon, explaining 17% of observed variability (Table 5(b)). In all cases, the topographic variables were better correlated with spatially smoothed predictions than those made from unsmoothed imagery. The comparatively poor correlation with topographic indices observed in this study might indicate that soil carbon and elemental nutrient content is more related to variability in field management and manure application than to elevation and soil redistribution, or that a more complex set of topography-influenced ecological and physical processes than was measured is at play. On low relief (<2% slope), ditch drained, Coastal Plain soils such as were sampled in this study, soil redistribution due to erosion is likely limited to short distances (<1–10 m), and prior converted wetland areas within fields can exhibit retained capacity for redoxomorphic activities that affect iron and aluminum transformations and increase denitrification and carbon accumulation. In this landscape, spatial distribution of soil carbon and nutrient content is likely influenced by a complex interaction of drainage status and land use history that is not easily characterized [5], and ultimately it may be more feasible to map surficial soil parameters using remote sensing technology than to predict them based on measurement and understanding of site specific processes.Table 5
Correlation between predicted analyte concentrations and topographic indices including (a) relative elevation data derived from 3 m LIDAR DEM, and (b) wetness index calculated from 3 m LIDAR DEM, for one of the six sampled fields (Mason). Predicted values were calculated for unsmoothed imagery (2.5 m pixel size) and for smoothed imagery to which a3×3 low band filter had been applied. See Section 2 for analyte descriptions.
CSiltAlFeR2R2R2R2(a) correlation with normalized elevationunsmoothed0.1000.0000.2320.013smoothed0.0510.0070.3860.000(b) correlation with wetness indexunsmoothed0.1830.0080.1610.091smoothed0.1780.0210.3170.037Soil moisture can become an interfering variable in remote sensing chemometric analyses, as it mutes the soil’s reflectance signal in a number of specific wavebands that are associated with water absorption bands [31]. Soil moisture for the 315 samples included in this study ranged between 5% and 42%, with most samples falling near the mean of 15% moisture.   When  observed moisture content was incorporated with the first derivative spectra of soil samples as a 177th predictor variable (following appropriate mean centering and variation normalization), the prediction accuracies associated with PLS models for C, silt, Al, and Fe did not improve (change in R2 of <0.02).  Similarly, when silt was included as an additional predictor variable, to possibly account for interferences between soil particle size distribution and spectral reflectance, the prediction accuracies associated with PLS models for C, silt, Al, and Fe did not improve. These results verify that PLS analysis is a robust method that can successfully extract predictive information from remotely sensed imagery of in situ agricultural soils.


## 3.1. Soil Test Results

Observed analyte concentrations for the six sampled fields are described in Table1. Soils were dry at the time of sampling, with moisture content ranging from 10% to 23% for the majority (91%) of samples (Table 1).  Overall, the distribution of observed carbon content values was not large, ranging from 0.6 to 2.0% (Table 1).  The Mason site, which had recently transitioned from long term conventional grain production to organic grain production, exhibited low soil carbon (mean 0.7%) content relative to the other fields (means of 1.2 to 1.4%).  The Schrader site had received long-term applications of dairy manure and exhibited somewhat increased soil C (mean 1.4%) relative to the other fields. The remaining four fields (Table 1: Temple 1S, Temple 1N, Temple 2, and Temple 3) were managed by one farmer, and were under similar management (conventional corn-wheat/soybean rotation). Overall, carbon content was low, as is typical of Maryland Eastern Shore farmland. Concentrations of the remaining analytes were within the normal range for agricultural soils, although variability among sites was not great.Table 1
Observed laboratory analyte concentrations for the six bare soil fields that were sampled on 10 April, 2007, with number of samples collected(n). See Section 2 for analyte descriptions.
FieldC%Sand%Silt%Clay%pHOM%K mg/kgCa mg/kgMg mg/kgMn mg/kgZn mg/kgFe mg/kgAl mg/kgN mg/kgP mg/kgCu mg/kgB mg/kgS mg/kgAl mg/kgMoisture%Elevation mWetnessTemplemin1.0736.528.93.34.61.70943785581.1973270.03231.80.2120.53485−5.07.21Smax1.6365.451.812.56.42.7097412142421116.94097980.125856.00.7982.0798352.712.7(n=46)Mean1.3248.642.49.05.72.22247817142342.52055010.07953.30.4041.4503150.010.27.3 hasd0.157.36.32.00.40.2515417139311.1751210.021100.80.1412.511941.51.1Templemin0.8932.324.65.84.71.30835016991.31152990.02221.90.1927.729911−4.87.11Nmax1.7566.955.716.66.33.0049312522571153.93609600.121755.00.7475.0960202.713.3(n=52)mean1.2149.540.99.65.32.12226785137322.52075590.06673.30.3746.3557160.010.47.1 hasd0.177.66.72.20.40.338919246290.7711460.02290.60.1413.914821.51.2Templemin0.8834.124.73.95.41.404348210740.81293060.03111.90.1726.630611−5.06.82max1.9168.255.417.37.03.105711423420224.44036480.112655.90.67108.8648421.512.7(n=55)mean1.2947.642.59.86.32.17164953256101.92314630.07493.30.4148.2463170.09.88.9 hasd0.227.66.22.40.40.351031856330.774680.02390.90.1319.27051.21.0Templemin1.0019.031.46.35.61.507161513481.4913520.03222.00.3523.83528−5.95.73max2.0459.567.919.77.23.4044518574881128.738010490.112439.14.42147.31049341.012.4(n=86)mean1.3232.954.912.26.32.50210992226343.01836520.07954.00.9337.1652140.010.418.1 hasd0.188.07.12.90.30.377027164231.2771880.01511.10.7714.518840.51.0Masonmin0.5547.69.51.35.70.504335775202.5715670.02342.20.2917.15676−2.76.4(n=30)max1.4087.743.410.96.82.202137752239412.213110330.0732021.20.6129.81033182.913.414.6 hamean0.7569.423.76.96.30.97111540133515.71028030.031484.40.4221.7803100.010.5sd0.199.77.92.50.30.394211134201.9171250.01713.40.083.412531.01.3Schradermin1.0937.218.06.04.92.1025347899224.11205410.04863.50.4428.15416−2.77.9(n=46)max1.8873.945.918.76.83.5060010662121118.64108220.0923212.90.8544.1822232.911.99.8 hamean1.4151.937.710.46.12.65449801160526.01776690.071336.60.6235.6669150.010.6sd0.188.66.72.70.40.348911327191.062660.01351.80.103.46631.00.5Overallmin0.5519.09.51.34.60.50433575540.8712990.02111.80.1717.12995−5.95.7(n=315)max2.0487.767.919.77.23.50974185748811512.241010490.1258521.24.42147.31049422.913.465.8 hamean1.2646.843.110.16.02.21236855186343.41905980.06934.10.5739.5598150.010.3sd0.2513.111.22.90.50.5613623870261.8771680.01531.30.3012.312741.11.0


## 3.2. Choice of Math Treatment

A SAS “shotgun approach” for hyperspectral data processing has previously been used to provide a factorial comparison of the effects of various math pretreatments on PLS analysis of soil spectra [23].  This method was used to compare a total of 15 different spectral pretreatments including 1st and 2nd derivatives with various gap widths (Table 2). Results showed that there were no statistical differences among the majority of the treatments, with the simplest treatments (no derivative, first derivative gap 2) often resulting in the best fit.  This led to the conclusion that the PLS data mining techniques are capable of extracting the majority of signal information from the untransformed spectral data, without the need for math pretreatment.Table 2
Partial least squares (PLS) prediction model goodness of fit (R2) associated with each of 15 math treatments, for the 13 analytes that predicted with R2>0.5, calculated using data from all 315 soil sampling locations1. The first derivative gap two (1Dg2, depicted in bold) was selected as the overall most successful model although it was occasionally outperformed (italic). See Section 2 for analyte descriptions.
DerivativeGapCSandSiltClaypHOMKCaMgMnZnFeAlaverageR2R2R2R2R2R2R2R2R2R2R2R2R2R2NON00.5780.7620.7630.5850.4420.6850.5550.6690.6920.6420.6700.7540.7770.6591ST10.5550.7700.7610.5650.5170.7060.5220.6540.7190.6240.6060.7070.7990.6541ST20.5910.7630.7630.6170.5490.7170.5140.6760.7080.6380.6470.7400.7820.6701ST40.5950.7540.7400.5960.4510.6920.5130.6300.6160.5650.6750.7030.7730.6391ST80.5840.7480.7390.6000.4130.6680.5780.6350.6420.5800.6360.6370.7270.6301ST160.5880.7400.7440.5500.4420.6720.5470.6100.6190.5670.5910.6660.7250.6201ST320.5250.7300.7220.5540.3810.6400.5200.6480.6190.5420.5800.6690.7650.6071ST640.3380.6230.6190.4700.3170.5370.4130.5190.5800.4880.5680.6360.6690.5212ND10.5420.5850.5970.4240.3610.5260.4840.7080.7050.5290.5990.6650.7290.5732ND20.5070.6810.6820.3720.4990.5660.4930.6650.6640.5660.5910.6730.7690.5952ND40.5760.7140.7240.4090.4490.6600.4840.6720.6690.6060.6130.7490.7520.6212ND80.5360.6980.6890.5330.4720.6020.5280.5950.6310.5420.6380.7430.7530.6122ND160.5800.7290.7250.5730.4870.6610.4980.5900.6400.5980.6170.6800.7120.6222ND320.5640.7150.6980.6090.4320.6310.5480.6220.5920.5480.6060.6390.6880.6072ND640.5250.5960.5890.4210.1560.5990.5110.5560.4640.5010.5380.6020.6140.5131Analytes that predicted poorly R2  (<0.5) included: N (<0.303), P (<.355), Cu (<0.358), B (<0.169), S (<0.282), and P saturation (<0.127).There is some argument to be made that different math treatments are more appropriate for particular analytes, due to the physical interaction of light with those particular constituents.  Some evidence for this was shown for potassium, for which a first derivative gap 8 increasedR2 from 0.514 (first derivative gap 2) to 0.578 (Table 2).  However, caution must be taken to prevent overfitting of datasets, and further study is needed to justify the selection of a diversity of math treatments for use with particular analytes.  Ultimately, the first derivative gap 2 was selected as the best math overall treatment even though it was slightly, but not significantly, outperformed in several cases (Table 2).


## 3.3. Partial Least Squares Regression: Calibration Dataset

Once the best math treatment was decided upon, PLS regression was performed in Matlab using a calibration dataset of 269 samples (all sampling locations from five of the six fields), with results reported in Table3. It should be noted that spectral data for the 315 samples used to calculate the SAS PLS in Table 2 were extracted from an earlier version of image processing output in which the 15 band errors had not yet been corrected and to which cross-track illumination had been applied. This discrepancy resulted in slightly lower R2 being used in Table 2 than in Table 3, but was not thought to have affected the relative performance of math treatments. Using Matlab PLS on 269 samples, 13 of the 19 analytes were predicted with R2>0.50 (Table 3(a)), and the remaining seven were predicted with R2<0.40: N (<0.30), P (<0.36), Cu (<0.36), B (<0.17), S (<0.28), and P saturation (<0.13).  A comparison of observed and predicted values for a selection of analytes (carbon, silt, aluminum, and iron) is shown in Figure 2.Table 3
Partial least squares (PLS) model accuracy in predicting soil analyte concentrations1 for (a) the 269 calibration samples using repeated leave-one-out cross validation and (b) 46 validation samples from the field (Temple 1S) that had been left apart from the calibration. Results were derived using a first derivative gap 2 math pretreatment, using spectra derived from unsmoothed imagery (2.5 m2 pixel size) or from spatially smoothed imagery to which a 3×3 low pass filter had been applied. See Section 2 for analyte descriptions. Units refer to analyte residual mean squared error (RMSE) values. Bold indicates sd/sepred >1.20.
Factors no.C%Sand%Silt%Clay%pHOM%K mg/kgCa mg/kgMg mg/kgMn mg/kgZn mg/kgFe mg/kgAl mg/kg8910910104891081010(a) 269 calibration samples collected from five agricultural fieldsusing unsmoothed imageryR20.650.790.790.660.510.750.590.690.690.620.640.750.76RMSE0.197.96.92.20.40.489.5166.150.319.61.449.3104.7using smoothed imageryR20.640.800.800.680.580.770.610.710.710.670.670.780.81RMSE0.187.26.22.00.40.386.1151.845.517.81.243.889.0(b) 46 validation samples collected from the remaining agricultural field (Temple 1S)using unsmoothed imagery2se pred0.209.07.32.40.50.4147.1178.151.427.51.252.9123.93bias0.10−2.881.421.000.030.3424.9240.4531.843.230.452.12107.33bias/mean (pred)0.08−0.060.030.110.010.150.100.050.220.090.180.010.214sd/se pred0.740.810.860.840.890.631.050.960.761.150.921.410.98using smoothed imageryse pred0.159.57.92.90.40.3154.9233.658.538.31.398.0170.4bias−0.10−0.87−0.660.95−0.25−0.13−23.06−2.46−2.07−0.380.10−5.2677.66bias/mean (pred)−0.08−0.02−0.020.10−0.05−0.06−0.100.00−0.02−0.010.04−0.030.14sd/se pred0.861.291.261.070.970.741.041.341.411.250.981.911.291Analytes that predicted poorly (R2 < 0.5) in the set of 315 samples (N, P, Cu,  B, S, and P saturation) are not included here2se pred: sd(observed-predicted)3bias: mean(pred)-mean(obs)4sd: standard deviation.Comparison of observed analyte concentrations with values predicted using partial least squares (PLS) regression on reflectance spectra extracted from smoothed imagery to which a3×3 low-pass filter had been applied, for (a) soil carbon, (b) silt, (c) iron, and (d) aluminum content.  The field containing 41 sampling points that were left out for validation (Temple 1S) is depicted with solid circles.
(a)(b)(c)(d)Prediction accuracy for carbon (R2=0.65) fell well within the range of results (0.45 to 0.98) found in a survey of 44 studies [20] but was somewhat poor in comparison with results found in some other studies for example, [11, 14, 15].  The somewhat poor R2 for carbon might lead to the conclusion that the sensor did not capture a good reflectance signal due to the increased noise often associated with airborne sensors that derives from atmospheric effects and variations in sensor-soil-sun geometry across the imagery. Indeed, a number of studies have found remote sensing spectroscopy to have reduced signal : noise relative to laboratory-based measurements [24].  However, the carbon PLS results reported in McCarty et al. [28], calculated for all 315 soil samples in this data set, using repeated leave-one-out correlation to determine goodness of fit, showed that spectra from the airborne sensor were as effective (R2=0.67) as spectra from a near-infrared benchtop laboratory spectrophotometer (R2=0.64) in predicting soil carbon concentrations.  For the remaining 12 analytes under consideration, the airborne sensor exhibited decreased accuracy relative to the laboratory spectrophotometer in three cases (change in R2 of −0.03 to −0.08) and exhibited increased accuracy in nine cases (change in R2 of 0.01 to 0.20).These observations led to the conclusion that the airborne sensor provides a viable option for mapping soil properties and that the somewhat poor prediction accuracies observed in this experiment apparently stemmed from features associated with the local soil environment or the calibration dataset, rather than the effectiveness of the airborne sensor. The distribution and range of observed analyte concentrations within a calibration dataset can have a substantial impact on prediction accuracy, and the somewhat low prediction accuracy for carbon content observed in this study (R2=0.65) is likely associated with the limited range and low values of soil carbon contents found within the calibration dataset (0.6% to 2.0%).For the other analytes, prediction accuracies (Table3(a)) were comparable with those found in other studies [8, 15, 17, 18]. Using principal components regression, Chang et al. [18] predicted Melich III analyte concentrations, and, similarly to this study, found Cu and P to have poor predictions, silt and clay content to predict with R2~0.8, and a number of other analytes, including Fe, Mg, Mn, K, and pH, predicting with R2>0.6.  It is not always known whether successful predictions are the direct action of the analyte upon the reflectance signal or instead the covariation of analyte concentrations with other factors which are influencing the spectral response [29]. Phosphorus, for example, has no expected reflectance resonance and is generally poorly predicted (e.g., 0.36 in this study, 0.40 in [18], and as low as 0.10 in other studies [8, 29]), and yet occasionally will be predicted with considerable accuracy [8, 29], likely due to covariation with spectrally responsive factors associated with labile organic matter.


## 3.4. Validation

Chemometric predictions of soil properties are typically validated by predicting analyte concentrations for samples that were not included in the calibration data set, and comparing the predicted results to observed values.  This can be achieved by calculating repeated leave-one-out or leave-ten-out predictions, or by selecting a distinct set of samples either randomly, spatially, or based on even sample distribution within the range of observed values [24, Section 18].  We chose to adopt a rigorous validation by removing an entire agricultural field comprising 46 samples (Temple 1S) from the data set, leaving 269 samples for calibration.  This field was one of four Temple fields that received similar crop management (corn-soy/wheat rotation with full tillage and no use of cover crops) and it exhibited analyte concentrations in the mid-range of the six fields (Table 1, Figure 2). When predicted values were calculated for the 46 validation samples and compared to observed values, bias ranged from 1% to 22% of observed analyte concentrations (Table 3(b)).  Seven of the analytes exhibited sd/se pred >1.20 (bold text in Table 3(b)) using the spectra extracted from smoothed imagery, indicating that the prediction was useful for those elements (sand, silt, Ca, Mg, Mn, Fe, and Al), while only one analyte (Fe) exhibited sd/se pred > 1.20  using the spectra extracted from unsmoothed imagery, indicating that spatial smoothing helped to increase signal : noise in the spectral data set.


## 3.5. Prediction Maps

For each analyte, a PLS regression vector of 178 coefficients was exported from Matlab and applied to the mean-centered first derivative reflectance spectra associated with each pixel of the hyperspectral imagery to calculate a 2.5-m raster map of predicted analyte concentrations. Prediction maps for a selection of analytes (carbon, silt, iron, and aluminum) were derived in this manner using both unsmoothed and smoothed imagery (Figure3). Ben-Dor et al. [19] have discussed the need for attention to imagery processing methods to extract maximum information from spectroscopic imagery. In this study, smoothing the hyperspectral imagery with a 3×3 low band filter prior to spectral data extraction resulted in somewhat improved prediction accuracies in comparison to unsmoothed imagery (Table 3(a)).  While this result implies that the smoothed extent of nine adjacent pixels (56.2 m2) provided a better average representation of soil characteristics at each sampling point than did individual 6.25 m2 pixels, it is more likely that the increased accuracy associated with the smoothed imagery is attributable to improvement of signal and reduction in sensor detector element variability (noise) that results from coaddition of the nine adjacent spectra.  Smoothing also generally improved the range and distribution of predicted analyte values found within each field image (Table 4), likely due to the correction of aberrations within the detector elements. However, it should be noted that discontinuities (striping) along the plant’s flight path were evident in both smoothed and unsmoothed predicted imagery (Figure 3), indicating the effects of variable detector element sensitivity along the sensor array and the need for continued improvement in sensor technology.  Several additional imagery adjustments were considered, including cross-track illumination, transformation from reflectance (R) to log(1/R), and also spectral smoothing (not shown). However, none of these operations improved results and they were therefore avoided.Table 4
Predicted concentrations of select analytes (C, Silt, Al, Fe)  derived fromn pixels of near infrared hyperspectral imagery for each of the six bare soil fields using either unsmoothed (1-pixel) imagery with 2.5 m resolution, or spatially smoothed (9-pixel) imagery to which a 3×3 low pass filter had been applied1. See Section 2 for analyte descriptions.
CSiltAlFeField1-pixel%9-pixel%1-pixel%9-pixel%1-pixel mg/kg9-pixel mg/kg1-pixel mg/kg9-pixel mg/kgTemple 1Smin0.971.1314.614.4342454056n=11731max1.861.9480.158.71017947474394mean1.411.4842.141.1633631194192sd0.100.105.95.098924644Temple 1Nmin0.850.961.05.6−311473764n=11436max1.952.0468.959.11304897452429mean1.231.2543.642.1589555200210sd0.110.115.84.91151126065Temple 2min−0.480.8114.68.734−3894054n=14221max2.072.9280.175.710171013474529mean1.271.2842.138.2633512194226sd0.250.185.97.1981464663Temple 3min0.390.600.21.2−1370−288−6n=28955max2.292.4880.874.911731031465474mean1.271.2752.351.7632639173170sd0.170.167.97.51541457163Masonmin0.410.44−7.6−7.5−19829000n=23351max2.271.9464.948.512311180391391mean0.810.8127.728.7755776111108sd0.110.106.47.1961003227Schradermin0.910.9812.313.42552965866n=15684max2.082.0561.656.9996971346315mean1.451.4838.638.3669681177172sd0.130.156.05.887813534Overallmin0.971.1314.614.42552965866n=105378max1.861.9461.648.5996897346315mean1.201.2241.340.56606481681711See Table 3 for the accuracy (R2) associated with each analyte prediction model.Maps of predicted soil carbon content calculated from (a) unsmoothed imagery (1-pixel data extraction) and (b) spatially smoothed imagery (9-pixel data extraction).  Two of the six sampled fields are depicted here, with Temple 3 on the left and Schrader on the right. Predicted values ranged from 0.4% (black) to 2.5% (white).
(a)(b)


## 3.6. Topographic Analysis

Visual interpretation of field topographic features relative to the predicted analyte maps revealed accumulation of carbon, silt, iron, and aluminum in low areas (Figure4).  This result is in agreement with previous studies in Iowa that have linked increased soil C carbon contents in low areas with soil redistribution processes [3], and with a study by Terra et al. [1] who detected a correlation of soil properties with topographic indices and electrical conductivity that explained 50% of observed variability in carbon content.  Other studies have evaluated the links between topography and soil constituents such as organic carbon at the landscape scale, showing accumulation at the bottom of slopes and in valley bottoms [30].  This study’s work was conducted at the field scale, reflecting microtopographic soil distribution processes within a relatively flat elevation gradient (0%–5% slope), as well as the formative effects of soil moisture distribution and spatially variable biological processes resulting from soil heterogeneity and the balance of redoxomorphic status between oxic and anoxic conditions. Understanding distribution patterns of soil analytes, and the influence of hydrogeomorphic controls on these patterns, can provide greater certainty about the influence of soil erosion and soil ecology on the fate of carbon and nutrients.Map of predicted values for select analytes (C, Silt, Fe, Al), overlaid on a high-resolution digital elevation map, for one of the six sampled fields (Mason). Predicted values were derived from imagery that had been smoothed with a3×3 low pass filter.
(a)(b)(c)(d)To investigate potential causes of spatial distribution in predicted analyte values, a 2.5 m elevation map (re-sampled from 3-m LIDAR DEM) was used to calculate a normalized elevation value (observed—field mean elevation) and a soil wetness index (two consecutive iterations of a SAGA enhanced lee3×3 filter) for each pixel of the predicted analyte maps. The normalized elevation variable was somewhat correlated (Table 5(a)) with predicted concentrations of aluminum (R2=0.23, 0.39 for unsmoothed and spatially smoothed imagery, respectively) and was not well correlated with predicted concentrations of silt, iron, and carbon.  The wetness index was similarly correlated with predicted aluminum concentration (R2=0.16,0.32) and was a better predictor for soil carbon, explaining 17% of observed variability (Table 5(b)). In all cases, the topographic variables were better correlated with spatially smoothed predictions than those made from unsmoothed imagery. The comparatively poor correlation with topographic indices observed in this study might indicate that soil carbon and elemental nutrient content is more related to variability in field management and manure application than to elevation and soil redistribution, or that a more complex set of topography-influenced ecological and physical processes than was measured is at play. On low relief (<2% slope), ditch drained, Coastal Plain soils such as were sampled in this study, soil redistribution due to erosion is likely limited to short distances (<1–10 m), and prior converted wetland areas within fields can exhibit retained capacity for redoxomorphic activities that affect iron and aluminum transformations and increase denitrification and carbon accumulation. In this landscape, spatial distribution of soil carbon and nutrient content is likely influenced by a complex interaction of drainage status and land use history that is not easily characterized [5], and ultimately it may be more feasible to map surficial soil parameters using remote sensing technology than to predict them based on measurement and understanding of site specific processes.Table 5
Correlation between predicted analyte concentrations and topographic indices including (a) relative elevation data derived from 3 m LIDAR DEM, and (b) wetness index calculated from 3 m LIDAR DEM, for one of the six sampled fields (Mason). Predicted values were calculated for unsmoothed imagery (2.5 m pixel size) and for smoothed imagery to which a3×3 low band filter had been applied. See Section 2 for analyte descriptions.
CSiltAlFeR2R2R2R2(a) correlation with normalized elevationunsmoothed0.1000.0000.2320.013smoothed0.0510.0070.3860.000(b) correlation with wetness indexunsmoothed0.1830.0080.1610.091smoothed0.1780.0210.3170.037Soil moisture can become an interfering variable in remote sensing chemometric analyses, as it mutes the soil’s reflectance signal in a number of specific wavebands that are associated with water absorption bands [31]. Soil moisture for the 315 samples included in this study ranged between 5% and 42%, with most samples falling near the mean of 15% moisture.   When  observed moisture content was incorporated with the first derivative spectra of soil samples as a 177th predictor variable (following appropriate mean centering and variation normalization), the prediction accuracies associated with PLS models for C, silt, Al, and Fe did not improve (change in R2 of <0.02).  Similarly, when silt was included as an additional predictor variable, to possibly account for interferences between soil particle size distribution and spectral reflectance, the prediction accuracies associated with PLS models for C, silt, Al, and Fe did not improve. These results verify that PLS analysis is a robust method that can successfully extract predictive information from remotely sensed imagery of in situ agricultural soils.


## 4. Conclusions

Aircraft-based acquisition of hyperspectral reflectance imagery can currently provide the necessary data to map soil properties in an efficient and rapid manner, and future improvements in sensor technology are expected to improve the signal : noise and spatial resolution associated with remote sensing imagery.  Results of this remote sensing study, conducted on relatively flat, moderately well drained, Coastal Plain silt loam soils, showed that the PLS framework is robust, and spectral math pretreatments can be kept simple, with a first derivative gap 2 providing good results for all analytes.  Spatial smoothing of reflectance imagery using a3×3 low-band filter improved results, likely by reducing noise inherent to variability within the spectrophotometer’s detector array. Fifteen of the 19 tested analytes predicted with R2>0.50 (with R2 from 0.51 to 0.79 for C, sand, silt, clay, pH, OM, K, Ca, Mg, Mn, Zn, Fe, and Al content; Table 3).  Prediction accuracy for carbon was not particularly good (R2=0.64) for this selection of sampled fields, perhaps owing to the low concentrations and limited range of observed variability (0.5% to 2%) in the calibration data set, but prediction accuracies for other analytes were on par with those found in other studies [8, 15, 18]. Predicted map values of select analytes, particularly aluminum, were correlated to field topography, indicating the influence of environmental processes on soil properties.  Further research is needed to understand the interplay of measurable topographic and hydrological variation inherent in this relatively flat agricultural landscape with management variables such as tillage, crop rotation, and additions of manures and fertilizers that also affect soil chemical and physical properties.Remote sensing approaches that use multispectral or hyperspectral imagery to predict soil properties are necessarily limited to analysis of plowed fields with very little crop residue or vegetation cover (bare soil). As a result, techniques that extrapolate soil properties data based on correlations derived from more easily measured parameters, such as topography, wetness, and cropping pattern may provide.  The fields used in this study were carefully selected to have negligible amounts of surface residue, green vegetation, or hydrologically active areas.  If this technique were put in practice on a larger scale, band filtering could be employed to remove such non-soil areas from analysis based on thresholding of spectral indices that correlate well with biomass (e.g., normalized difference vegetation index [32]), crop residue (e.g., cellulose adsorption index [33]), or soil moisture (e.g., water adsorption [13, 31]).As hyperspectral imagery becomes more readily available and at lower cost, the application of partial least squares (PLS) regression to soil spectral reflectance data can provide an effective method for calculating high-resolution raster maps of important soil properties including texture, pH, and carbon and nutrient content.  This information can then be used to inform farmer decision making, support precision environmental management of agricultural fields, increase sustainable crop production, and help to reduce nutrient, sediment, and carbon losses from agricultural systems.



---

*Source: 358193-2011-08-11.xml*