Linear Regression Trend Analysis
To perform a simple Linear Regression, we first need to access and prepare the Image and Geometry data for our Area of Interest.
In this example, we are going to use all available Landsat 8 Tier 1 Surface Reflectance data with a cloud coverage of less than 10% for the region of Homs, Syria.
Preparing the data
//Access the Area of Interest geometry
var extent_aoi = ee.FeatureCollection("FAO/GAUL/2015/level1")
.filter(ee.Filter.eq('ADM1_NAME', 'Homs')); //filter for entry that equals the ADM1_NAME 'Homs'
print(extent_aoi, 'Extent AOI');
Map.addLayer(extent_aoi, {},'Extent AOI');
Map.centerObject(extent_aoi, 8);
//Access Landsat 8 data from the day it launched until today. Filter for the desired region and less than 10% cloud cover.
var ls8 = ee.ImageCollection('LANDSAT/LC08/C01/T1_SR')
.filterBounds(extent_aoi)
.filter(ee.Filter.lt('CLOUD_COVER', 10));
print(ls8, 'Image Collection 2013-2021 <10% clouds'); //Use the console to look at the amount of Images accumulated
var visParams = {'min': 400,'max': [4000,3000,3000], 'bands':'B5,B4,B3'}; //bands 5, 4 and 3 are to be visualized; False Color Composite
Map.addLayer(ls8, visParams ,'Image Collection 2013-2021 <10% clouds', false); // Visualize the ImageCollection in the Map panel and assign a name to it; Hide the layer by default, as it takes long to render
Add Time- and NDVI-Bands
Next, we need to create our dependent and independent variables. In this example, we want to analyse the trend of NDVI-Values over time. This means that we need to write two functions that allow us to map the NDVI and time information to every Image of our of our ImageCollection. We can use .rename to name the added NDVI-band 'NDVI' , which makes browsing the data easier .
var addNDVI = function(image) {
var ndvi = image.normalizedDifference(['B5', 'B4']).rename('NDVI');
return image.addBands(ndvi);
};
//Apply the function to the ImageCollection
var ls8_ndvi = ls8.map(addNDVI);
print(ls8_ndvi, 'Image Collection with NDVI')
//Step 3: Write a function that adds a time band to the image collection. The time is divided by a large number to avoid very small slopes in the regression output, hence making it easier to visualize and interpret.
var createTimeBand = function(image) {
return image.addBands(image.metadata('system:time_start').divide(1e15));
};
//Apply the function to the ImageCollection.
var ls8_trend = ls8_ndvi.map(createTimeBand);
//Have a look at the prepared ImageCollection
print(ls8_trend, 'Image Collection with NDVI and time band');
Run a Linear Fit Function
Now that we have our image data and variables ready, we can easily calculate the linear fit of the regression using the reducer .linearFit() to analyse the NDVI-trend of Homs since Landsat 8 started recording. Check its documentation at Docs -> ee.Reducer -> ee.Reducer.linearFit to find out more about its syntax. The first variable is the independent variable, in our case the time-band, and the second variable the dependent, in our case the NDVI-band.
//Calulate the NDVI-Trends using the linear fit reducer. The first variable is the independent variable, the second one the dependent.
var linearFit_ndvi = ls8_trend.select(['system:time_start', 'NDVI'])
.reduce(ee.Reducer.linearFit());
//Have a look at the output in the Console. It consists of two bands, the "scale", which is the slope, and the offset of the linear regression.
print(linearFit_ndvi, 'Linear Regression Fit')
//Display the results. In this example, we will display it in two different ways: One showing the trends in colour, and one showing the trends in a greyscale.
//For both variants, we need to use the Inspector to evaluate a suitable range of max-values to visualize the data.
//Trends displayed in colour
Map.addLayer(linearFit_ndvi.clip(extent_aoi),
{min: 0, max: [-130, 0.9, 130], bands: ['scale', 'offset', 'scale']}, 'Linear Fit NDVI colours');
//Negative slopes appear in red, showing negative trends -> vegetation loss!
//Positive slopes appear in blue, showing positive trends -> vegetation gain!
//Trends displayed in greyscale
Map.addLayer(linearFit_ndvi.select('scale').clip(extent_aoi),
{min:-130, max:130}, 'Linear Fit NDVI Greyscale');
//Negative slopes appear in black, showing negative trends -> vegetation loss!
//Positive slopes appear in white, showing positive trends -> vegetation gain!