Relationships that vary in space and time: A challenge for simple linear regression

Regression is a fundamental analytical tool for any researcher interested in relating cause and effect. A basic assumption of regression is that of stationarity, the principal that the relationship between a predictor variable and its response is constant across its sample space; that a relationship is true in all regions where it is being applied. This assumption is a particularly poor assumption in spatially based analyses, where we know that interactions may exist between known and often unknown factors in how response variable relates to a given explanatory variable. While this is a  challenge to simple linear regression, it is also what generally makes spatial problems interesting: the fact that relationships are not constant across space and time.

Spatially weighted regression challenges the assumption of stationary in that where simple linear regression develops a single relationship to describes a phenomena, spatially weighted regression allows the relationship to vary spatially. Unhinging the relationship between a explanatory variable and its response spatially creates a set of local coefficient for each instance where an explanatory variable is offered. This is done through the use of a  weighting function. Wherein simple linear regression, each data point assumes equal weight with regards to the final relationship, a weighting function applies greater import to values closer to where a regression point would be calculated.

---

Fig 1: A spatial weighting function weights data points closer to a  regression point. In this way bandwidths can vary across a feature space, such that two local regression values may be constructed of a different number of data points.

---

Fig 2: Spatially weighted regression allows the relationship between a response and explanatory to vary across a study region.

 

NDVI and weed density: A case study in spatially weighted regression

Normalized difference vegetation index (NDVI) has been used in remote sensing  as a proxy for phenology in many remote sensing and cropping systems studies. NDVI is calculated as the ratio of red to near-infrared light, and is generally related to the amount of green photo-synthetically active tissue. In principal, weeds and crops should be distinguishable based on how their NDVI response varies in time.

---

Question: Can NDVI be used to predict weed density? Does the relationship between NDVI and weed density vary spatially? 

---

Fig 3: A general hypothesis for how weeds and crop may in their NDVI response. Weeds may mature earlier or later than a crop, but this relationship may also vary spatially.

Here, spatially weighted regression is offered as a method for distinguishing relationships between weed density and NDVI. Allowing the relationship between NDVI and weed density to vary spatially may allow one classify areas of the field based on the nature of these relationships.

---

Fig 4: NDVI over the course of the growing season at a 17 acre field located NE of Pendleton, OR, Summer 2015. Note that the field does not increase or decrease in NDVI evenly, but rather peak NDVI passes as a wave across the field. [ Note: GIF does not appear animated in preview. Click on the GIF directly if animation is not working ]

An important first step in deciding if a data set is suitable for spatially weighted regression is to look at the residuals of the linear relationship you are choosing to model. Here we examine the following function for predicting weed density:

This function uses 3 samples of NDVI over the course of the growing season, centered around the time of peak NDVI. The purpose of using these 3 times was to try and emphasize time periods where weeds and crop would vary in their relative response in NDVI.

---

Fig 5: Weed densities were calculated based on linear transects made prior to harvest. Weed hot spots were calculated using the Getis-Ord Gi* statistic in ArcMap. Weed hot spots from a prior analysis were used as input for predicting weed density in this exercise.

---

Fig 6: Predicted weed density for a multiple regression model based on 3 measurements of NDVI surrounding peak NDVI.  Multiple R2:  0.044, P-value < 0.001

 

---

Fig 7: Residuals from a multiple regression model based on 3 measurements of NDVI surrounding peak NDVI.  Both high and low residuals cluster spatially, indicating that the relationship between NDVI and weed density may vary spatially and may be a good candidate for geographically weighted regression.

---

Fig 8: Predicted weed density from a spatially weighted regression. Quasi-global R2: 0.92.

By using a spatially weighted regression, we’re able to account for 92 percent of the variance that occurs in the distribution of weeds in this field. Unlike in a standard regression, the result of this process is a collection of local regression formula. In this sense, the result is not a result that can be easily extrapolated to predict weeds distributions in future data sets. However, these coefficients do offer us the opportunity to look for some spatial patterns that may yield additional information as to what the nature of these local spatial relationships might be

---

 

Fig 10:  Map classified by coefficient slope.

Introduction:

Hotspot analysis is a method of statistical analysis for spatially distributed data whereby both the magnitude of a point and its relationship to the points around it are taken into consideration to identify statistically significant extreme values in a data set. While a single spatial data point may represent an extreme value, in the context of weeds management what often matters is being able to identify aggregates of values which taken as a whole represent an extreme condition. Where one large weed may have some impact, many medium and smaller weeds are likely to have a disproportionate impact. This has implications for weed management and weed science, specifically of the variety of techniques available to users for mapping the spatial distribution of invasive weeds, most rely on expert knowledge to decide what constitutes a “weedy” versus a “non-weedy” condition. This consideration is made independently of the geographic context of a value, usually by visually evaluation of the distribution of measured values, and deciding some cutoff, real or perceived, in the distribution of values. These cutoffs are often subjective, and are typically yield poor results when applied to independently generated data sets. Hotspot  analysis offers us a road forward where statistically significant extreme regions can be identified, and various mapping techniques compared by how well they predict these extreme regions.

My goal for this project was to see if hotspot analysis provided a better technique for mapping weeds populations in wheat fields at the time of harvest. Previously, we visually evaluated the histograms, and using some knowledge of how the weeds were distributed to determine cutoffs for “weedy” and “non-weedy” values. Here I wanted to see if I could generate more accurate maps of weed density relative to a reference data set.

Sources of data:

Fig. 1: Sampling for weed density was done at the time of harvest using a spectrometer in the grain stream of a combine, and with a visual evaluator watching for green weeds during harvest.
In Summer 2015 we harvested a field winter wheat NE of Pendleton, OR using combine outfitted with a spectrometer in the grain. A ratio of red to near infrared reflectance was taken, along with visual observations made by an observer from inside the cab of the combine. Prior to harvest, we walked the field on gridded transects, recording observations of all weeds at the species level. These data were all re-gridded to a common 7m X 7m grid using bilinear interpolation, and brought into ArcGIS for analysis.

Classification:

In ArcGIS, hotspot analysis was conducted using Hotspot analysis tool. The Getis-Ord Gi* statistic identifies statistically significant hot and cold spots in whatever spatially referenced value you use as an feature class for classification. For the purposes of this exercise, positive Getis-Ord Gi* statistics were considered to be “weedy”, while not-significant and negative scores were considered to be “non-weedy”. Cutoffs based on a visual evaluation of histogram, and knowledge of what should represent a weedy condition used to distinguish between “weedy” and “non-weedy” values  and their associated maps classified accordingly.

 

Fig 2: Ground reference data set classified using Hotspot analysis and histogram classification. All green values were considered positive for weeds.

Fig 3: Spectral assessment of weeds distribution classified using Hotspot analysis and histogram classification. All green values were considered positive for weeds.

 

Fig 4: Visual assessment of weeds distribution classified using Hotspot analysis and histogram classification.  All green values were considered positive for weeds.

Comparing techniques
Classified maps were then compared for accuracy using a confusion matrix approach. Accuracy here is defined as the rate of true positives and true negatives, divided by the total number of samples. Classified maps were also compared for their precision and sensitivity, where precision is the rate of true positives of the condition positives, and precision is the rate of true positives over the predicted positives. Precision addresses the question of if I predict it positive, how often am I correct? Sensitivity addresses the question, when the reference data predicts it positive, how often was my prediction correct?

Fig 5: Error assessment comparing spectral assessment of weediness to ground reference data using the histogram classification.

Fig 6: Results of the error assessments for all mapping techniques compared with their ground reference data.

Results:

These results show that hotspot analysis was a superior method for image classification when compared with a visual evaluation of the distribution.  This is probably a result of the fact that hotspot analysis takes into consideration the neighborhood of values a measurement resides in. Hotspot analysis increased the accuracy for both the spectral and the on-combine visual assessments. There were mixed results regarding its impacts on sensitivity and precision however. Hotspot classification increased the precision of the visual assessment, but did so at a cost to sensitivity. Conversely, the hotspot technique increased the sensitivity of the spectral assessment, but did so at a cost to its precision. It may be that precision comes at a cost to sensitivity using this method for comparing classifications. In general however, the most substantial gains were in accuracy, which is the most important factor for this mapping effort.

Question:

Foundational to effectiveness of management of invasive species in agricultural and native landscapes is the question of spatial extent and the ability to quantify the impacts of invasive species on agronomic efforts and native vegetation. Invasive plants represent a threat to both agricultural and native landscapes in the form of reduced ecosystem function, increased resource consumption, and reduced yields from agricultural systems. In spite of significant efforts to control and reduce impacts from invasive species, invasive weeds cause an estimated loss of $2 billion annually in the US. Currently, interspecific competition is one of the major limitations for oilseed  and grain production in dryland cropping systems of the Pacific Northwest (PNW). Opportunities for the for the precision monitoring of managed and native ecosystems have become available through the use of low altitude remote sensing systems and high resolution satellite systems. However, methods for resolving species level classification in high resolution multispectral remote sensing systems remain lacking. This is partially due to the relative novelty of these systems, but is also related the lack of suitable reference data at spatial and temporal scales for regionally based models. The broad research question I’m is how does spectral trajectory relate to weed density, and can this information be used to distinguish the spatial extent of weeds in dryland cropping systems? My prediction, is that by increasing the spatial and temporal resolution of these data, crop and non-crop species will be distinguishable based on their relative rate of change in greenness.

The objectives I have for this class are to 1) determine the spatial and temporal resolutions at which weed species are distinguishable from crop species using a spectral trajectory technique, 2) compare these methods with ground reference data in a dryland cropping systems study. The major outcome of this work would be a method for distinguishing weed species from crop species in a dryland environment, and the identification of the minimum temporal resolution for distinguishing species in multispectral imagery.

Dataset:

The data set I have for addressing these questions is a composite of 7 flights of images taken with a multispectral camera in a cropping systems study in Eastern Oregon. Flights were conducted in conjunction with visual estimates of weed density in semi-permanent monitoring frames installed into the cropping systems study. The images are currently at a low level of processing. One of my goals as a part of this project will be to orthomosaic the images such that I can perform a time series analysis across image collection dates. The temporal resolution is from 3-20 days between flights, whiles the spatial resolution is 3 cm. The images cover the entire spatial extent of the experiment.

Hypothesis:

What I plan to do for this experiment is that after generating an orthoraster, I will be able to distinguish between the quantify of weed species and crop species in a frame based on the spectral trajectory of individual pixels in that frame. The question I hope to answer will be how distant in time do sample dates have to be before weed species can be distinguished from crop species based on their spectral trajectory.

Approach:

I plan on using a trial version of Agisoft to orthomosaic my images. It may be possible to conduct this analysis without performing an ortho mosaicing of the images, however, an orthomosaic will have a number of advantages to non-mosaiced images. Without an orthomosaic every individual analysis will have to be hard coded, whereas with an orthomosaic, I can automate much of the processing of these data.

Outcome:

The goal of this analysis will be to identify the minimum time required to discern species, and to describe the statistical relationship between species abundance based on spectral trajectory and ground reference data.

Significance:

While there has been a significant surge in the interest of UAV’s and low altitude remote sensing, the actual number of useful products for land managers to make decisions based on these data is very minimal. This work would identify the minimum temporal resolution a resource manager would need to have before they can identify weed species based on spectral characteristics.

Preparation:

I would say I have minimal experience in Arc-info, model builder and Python. I have moderate to high level experience in R. I would also consider myself to have moderate to high levels of experience in image processing and image analysis.