0

GEE, Sentinel and linear Regression

I am using GEE (Google Earth Engine) to compute a Sentinel-2 based change analysis. What I would like is basically a linear regression per pixel, with the slope as an output band. However, I need to address issues with atmospheric contamination of imagery. Hence, I am running an image classification to then produce a mask based on my classes for atmospheric contamination (not happy with QA60). I then tried to use

a) the linearFit reducer

Problem: If data is masked in only one scene within the collection, the linearFit is masked too (even though there might be other imagery available without mask)! However, I would like the masked data simply to be ignored as Null and the trend to be displayed using the rest of the available data. Hence I tried b), where beta can be specified as NULL

b) the rubustlinearRegression

Problem: I do not understand the output Layer! what is the weird crosses? If I use the pixel inspector it shows the blank regions as masked. Why is it masked differently? Is it because of the min and max values I specified for the lrImage layer? If that is the case, Is there a way to see the range of values within a band as to then define min and max correctly?

Generally, I am assuming that I cannot map the trend with clouds as noData because the cloud mask is run through the entire collection, which includes the time band. Would it be a solution to exclude the time band from the mask? If yes, how would I do that?

This is the link: https://code.earthengine.google.com/1cbaf9ad6265c166bdbe508c62d7ad4b

This is the code. It includes a point geometry = roi5, Sentinel2 collection = MSI, NSW = Polygon used as boundary, and Feature Collections as landcover classes.

//Compute a trend for the cloudfreecollection

// Select the bands to model with the independent variable first.
var trend = QAclclear.select(['system:time_start', 'NDVI'])
  // Compute the linear trend over time.
  .reduce(ee.Reducer.linearFit());

//the max: [-1, 1, 10] impacts the masking!!!!
Map.addLayer(
    trend,
    {min:0 , max: [-1, 1, 10000] , bands: ['scale', 'scale', 'offset']},
    'NDVI trend');

//Robust Linear Regression where beta (here Null) is the error outlier, and a default

// This function adds a constant band to the image.
var createConstantBand = function(image) {
  return ee.Image(1).addBands(image);
};

//Map the new constant band over the QAclclear (atmosph. dist. masked) collection, and select the constant, independent, and dependent value
var ConstantQAclclear = QAclclear
              .map(createConstantBand)
              // Select the predictors and the responses.
  .select(['constant', 'system:time_start', 'NDVI']);


print(ConstantQAclclear);


//Take the Collection with the constant value and apply the reducer
var RLR = ConstantQAclclear.select(['system:time_start', 'NDVI', 'constant'])
  // Compute the linear trend over time.
  .reduce(ee.Reducer.robustLinearRegression({numX: 2, numY: 1}));    

print(RLR);

// The results are array images that must be flattened for display.
// These lists label the information along each axis of the arrays.
var bandNames = [['constant', 'time'], // 0-axis variation.
                 ['NDVI']]; // 1-axis variation.

// Flatten the array images to get multi-band images according to the labels.
var lrImage = RLR.select(['coefficients']).arrayFlatten(bandNames);

// Display the OLS results.
Map.addLayer(lrImage,
  {min: 0, max: [1, 1], bands: ['constant_NDVI', 'time_NDVI']}, 'OLS');

I am new to this.

  • I would recommend to join the GEE Google group. – Thomas Jan 2 '18 at 2:27
  • 1
    Is the entire classification section of your code unrelated your question? If so, please remove it. Also note that the output of a linear regression reducer is an array image as described here. You need to flatten it. This tutorial may be useful for understanding how to manipulate array images. – Nicholas Clinton Jan 2 '18 at 19:13
  • 1
    You have included six questions in your post. Best practices are to ask a single question in each post. See: gis.meta.stackexchange.com/questions/3349/… – Tyler Erickson Jan 2 '18 at 19:37
1

To answer your first question:

I do not understand the output Layer! what is the weird crosses?

The crosses are caused by how Sentinel-2 images are distributed. There is repeated information along the edge of the "granules". This can be easily seen by displaying a count of the pixels:

var s2_granules = ee.ImageCollection('COPERNICUS/S2')
                   .filterDate('2017-01-01', '2017-01-02');

Map.addLayer(
  s2_granules.count(),
  {bands:'B1', min:0, max:6, palette:"black,green,grey,yellow,blue,white,red"},
  's2_granules.count()'
);

Code link: https://code.earthengine.google.com/d18a618c9b529eac197f46624dd4ad86

enter image description here

As can been seen by the red areas, information can be repeated in as many as 6 Sentinel-2 granules.

For more information, see the "Granules and Tiles" section of the Sentinel-2 User Handbook.

| improve this answer | |
  • @Caroline if the answer helped and/or answered the single question you wanted to ask the most then please upvote and/or accept it. – PolyGeo Jan 3 '18 at 7:02

Not the answer you're looking for? Browse other questions tagged or ask your own question.