Regressione lineare

Earth Engine dispone di diversi metodi per eseguire la regressione lineare utilizzando i riduttori:

  • ee.Reducer.linearFit()
  • ee.Reducer.linearRegression()
  • ee.Reducer.robustLinearRegression()
  • ee.Reducer.ridgeRegression()

Il riduttore di regressione lineare più semplice è linearFit(), che calcola la stima dei minimi quadrati di una funzione lineare di una variabile con un termine costante. Per un approccio più flessibile alla creazione di modelli lineari, utilizza uno dei riduttore di regressione lineare che consentono un numero variabile di variabili indipendenti e dipendenti. linearRegression() implementa la regressione con minimi quadrati ordinari(OLS). robustLinearRegression() utilizza una funzione di costo basata su residui di regressione per attenuare in modo iterativo gli outlier nei dati (O'Leary, 1990). ridgeRegression() esegue la regressione lineare con regolarizzazione L2.

L'analisi di regressione con questi metodi è adatta per ridurre gli oggetti ee.ImageCollection, ee.Image, ee.FeatureCollection e ee.List. Gli esempi riportati di seguito mostrano un'applicazione per ciascun caso. Tieni presente che linearRegression(), robustLinearRegression() e ridgeRegression() hanno tutti le stesse strutture di input e output, ma linearFit() si aspetta un input a due bande (X seguito da Y) e ridgeRegression() ha un parametro aggiuntivo (lambda, facoltativo) e un output (pValue).

ee.ImageCollection

linearFit()

I dati devono essere configurati come immagine di input a due bande, dove la prima banda è la variabile indipendente e la seconda è la variabile dipendente. L'esempio seguente mostra la stima della tendenza lineare delle precipitazioni future (dopo il 2006 nei dati NEX-DCP30) previste dai modelli climatici. La variabile dipendente è la precipitazione prevista e la variabile indipendente è il tempo, aggiunto prima di chiamarelinearFit():

Editor di codice (JavaScript)

// This function adds a time band to the image.
var createTimeBand = function(image) {
  // Scale milliseconds by a large constant to avoid very small slopes
  // in the linear regression output.
  return image.addBands(image.metadata('system:time_start').divide(1e18));
};

// Load the input image collection: projected climate data.
var collection = ee.ImageCollection('NASA/NEX-DCP30_ENSEMBLE_STATS')
  .filter(ee.Filter.eq('scenario', 'rcp85'))
  .filterDate(ee.Date('2006-01-01'), ee.Date('2050-01-01'))
  // Map the time band function over the collection.
  .map(createTimeBand);

// Reduce the collection with the linear fit reducer.
// Independent variable are followed by dependent variables.
var linearFit = collection.select(['system:time_start', 'pr_mean'])
  .reduce(ee.Reducer.linearFit());

// Display the results.
Map.setCenter(-100.11, 40.38, 5);
Map.addLayer(linearFit,
  {min: 0, max: [-0.9, 8e-5, 1], bands: ['scale', 'offset', 'scale']}, 'fit');

Tieni presente che l'output contiene due bande, "offset" (intercetta) e "scala" ("scala" in questo contesto si riferisce alla pendenza della linea e non deve essere confusa con il parametro di scala inserito in molti riduttori, ovvero la scala spaziale). Il risultato, con aree di tendenza in aumento in blu, tendenza in calo in rosso e nessuna tendenza in verde, dovrebbe assomigliare alla Figura 1.


Figura 1. L'output di linearFit() applicato alle precipitazioni previste. Le aree in cui si prevede un aumento delle precipitazioni sono mostrate in blu, mentre quelle in cui si prevede una diminuzione sono mostrate in rosso.

linearRegression()

Ad esempio, supponiamo che esistano due variabili dipendenti: precipitazioni e temperatura massima, e due variabili indipendenti: una costante e il tempo. La raccolta è identica a quella dell'esempio precedente, ma la banda costante deve essere aggiunta manualmente prima della riduzione. Le prime due bande dell'input sono le variabili "X" (indipendenti) e le due bande successive sono le variabili "Y" (dipendenti). In questo esempio, ottieni prima i coefficienti di regressione, poi appiattisci l'immagine dell'array per estrarre le bande di interesse:

Editor di codice (JavaScript)

// This function adds a time band to the image.
var createTimeBand = function(image) {
  // Scale milliseconds by a large constant.
  return image.addBands(image.metadata('system:time_start').divide(1e18));
};

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

// Load the input image collection: projected climate data.
var collection = ee.ImageCollection('NASA/NEX-DCP30_ENSEMBLE_STATS')
  .filterDate(ee.Date('2006-01-01'), ee.Date('2099-01-01'))
  .filter(ee.Filter.eq('scenario', 'rcp85'))
  // Map the functions over the collection, to get constant and time bands.
  .map(createTimeBand)
  .map(createConstantBand)
  // Select the predictors and the responses.
  .select(['constant', 'system:time_start', 'pr_mean', 'tasmax_mean']);

// Compute ordinary least squares regression coefficients.
var linearRegression = collection.reduce(
  ee.Reducer.linearRegression({
    numX: 2,
    numY: 2
}));

// Compute robust linear regression coefficients.
var robustLinearRegression = collection.reduce(
  ee.Reducer.robustLinearRegression({
    numX: 2,
    numY: 2
}));

// 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.
                 ['precip', 'temp']]; // 1-axis variation.

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

// Display the OLS results.
Map.setCenter(-100.11, 40.38, 5);
Map.addLayer(lrImage,
  {min: 0, max: [-0.9, 8e-5, 1], bands: ['time_precip', 'constant_precip', 'time_precip']}, 'OLS');

// Compare the results at a specific point:
print('OLS estimates:', lrImage.reduceRegion({
  reducer: ee.Reducer.first(),
  geometry: ee.Geometry.Point([-96.0, 41.0]),
  scale: 1000
}));

print('Robust estimates:', rlrImage.reduceRegion({
  reducer: ee.Reducer.first(),
  geometry: ee.Geometry.Point([-96.0, 41.0]),
  scale: 1000
}));

Esamina i risultati per scoprire che l'output linearRegression() è equivalente ai coefficienti stimati dal riduttore linearFit(), anche se l'output linearRegression() ha anche coefficienti per l'altra variabile dipendente tasmax_mean. I coefficienti di regressione lineare robusti sono diversi dalle stime OLS. L'esempio confronta i coefficienti dei diversi metodi di regressione in un punto specifico.

ee.Image

Nel contesto di un oggetto ee.Image, i riduttori della regressione possono essere utilizzati con reduceRegion o reduceRegions per eseguire la regressione lineare sui pixel delle regioni. I seguenti esempi mostrano come calcolare i coefficienti di regressione tra le bande Landsat in un poligono arbitrario.

linearFit()

La sezione della guida che descrive i grafici dei dati degli array mostra un grafico a dispersione della correlazione tra le bande SWIR1 e SWIR2 di Landsat 8. Qui vengono calcolati i coefficienti di regressione lineare per questa relazione. Viene restituito un dizionario contenente le proprietà 'offset' (intercetta sull'asse y) e 'scale' (pendenza).

Editor di codice (JavaScript)

// Define a rectangle geometry around San Francisco.
var sanFrancisco = ee.Geometry.Rectangle([-122.45, 37.74, -122.4, 37.8]);

// Import a Landsat 8 TOA image for this region.
var img = ee.Image('LANDSAT/LC08/C02/T1_TOA/LC08_044034_20140318');

// Subset the SWIR1 and SWIR2 bands. In the regression reducer, independent
// variables come first followed by the dependent variables. In this case,
// B5 (SWIR1) is the independent variable and B6 (SWIR2) is the dependent
// variable.
var imgRegress = img.select(['B5', 'B6']);

// Calculate regression coefficients for the set of pixels intersecting the
// above defined region using reduceRegion with ee.Reducer.linearFit().
var linearFit = imgRegress.reduceRegion({
  reducer: ee.Reducer.linearFit(),
  geometry: sanFrancisco,
  scale: 30,
});

// Inspect the results.
print('OLS estimates:', linearFit);
print('y-intercept:', linearFit.get('offset'));
print('Slope:', linearFit.get('scale'));

linearRegression()

Qui viene applicata la stessa analisi della sezione linearFit precedente, solo che questa volta viene utilizzata la funzione ee.Reducer.linearRegression. Tieni presente che un'immagine di regressione è costituita da tre immagini separate: un'immagine costante e le immagini che rappresentano le bande SWIR1 e SWIR2 della stessa immagine Landsat 8. Tieni presente che puoi combinare qualsiasi insieme di bande per creare un'immagine di input per la riduzione delle regioni per ee.Reducer.linearRegression, che non devono necessariamente appartenere alla stessa immagine di origine.

Editor di codice (JavaScript)

// Define a rectangle geometry around San Francisco.
var sanFrancisco = ee.Geometry.Rectangle([-122.45, 37.74, -122.4, 37.8]);

// Import a Landsat 8 TOA image for this region.
var img = ee.Image('LANDSAT/LC08/C02/T1_TOA/LC08_044034_20140318');

// Create a new image that is the concatenation of three images: a constant,
// the SWIR1 band, and the SWIR2 band.
var constant = ee.Image(1);
var xVar = img.select('B5');
var yVar = img.select('B6');
var imgRegress = ee.Image.cat(constant, xVar, yVar);

// Calculate regression coefficients for the set of pixels intersecting the
// above defined region using reduceRegion. The numX parameter is set as 2
// because the constant and the SWIR1 bands are independent variables and they
// are the first two bands in the stack; numY is set as 1 because there is only
// one dependent variable (SWIR2) and it follows as band three in the stack.
var linearRegression = imgRegress.reduceRegion({
  reducer: ee.Reducer.linearRegression({
    numX: 2,
    numY: 1
  }),
  geometry: sanFrancisco,
  scale: 30,
});

// Convert the coefficients array to a list.
var coefList = ee.Array(linearRegression.get('coefficients')).toList();

// Extract the y-intercept and slope.
var b0 = ee.List(coefList.get(0)).get(0); // y-intercept
var b1 = ee.List(coefList.get(1)).get(0); // slope

// Extract the residuals.
var residuals = ee.Array(linearRegression.get('residuals')).toList().get(0);

// Inspect the results.
print('OLS estimates', linearRegression);
print('y-intercept:', b0);
print('Slope:', b1);
print('Residuals:', residuals);

Viene restituito un dizionario contenente le proprietà 'coefficients' e 'residuals'. La proprietà 'coefficients' è un array con dimensioni (numX, numY); ogni colonna contiene i coefficienti per la corrispondente variabile dipendente. In questo caso, l'array ha due righe e una colonna: la riga 1, colonna 1 è l'intercetta sull'asse y e la riga 2, colonna 1 è la pendenza. La proprietà 'residuals' è il vettore della radice quadrata media dei residui di ogni variabile dipendente. Estrai i coefficienti trasmettendo il risultato come array e poi estraendo gli elementi desiderati o convertendo l'array in un elenco e selezionando i coefficienti in base alla posizione dell'indice.

ee.FeatureCollection

Supponiamo che tu voglia conoscere la relazione lineare tra la riflettanza SWIR1 di Sentinel-2 e quella di Landsat 8. In questo esempio, per calcolare la relazione viene utilizzato un campione casuale di pixel formattato come raccolta di elementi di punti. Viene generato un grafico a dispersione delle coppie di pixel insieme alla linea di migliore approssimazione per i minimi quadrati (Figura 2).

Editor di codice (JavaScript)

// Import a Sentinel-2 TOA image.
var s2ImgSwir1 = ee.Image('COPERNICUS/S2/20191022T185429_20191022T185427_T10SEH');

// Import a Landsat 8 TOA image from 12 days earlier than the S2 image.
var l8ImgSwir1 = ee.Image('LANDSAT/LC08/C02/T1_TOA/LC08_044033_20191010');

// Get the intersection between the two images - the area of interest (aoi).
var aoi = s2ImgSwir1.geometry().intersection(l8ImgSwir1.geometry());

// Get a set of 1000 random points from within the aoi. A feature collection
// is returned.
var sample = ee.FeatureCollection.randomPoints({
  region: aoi,
  points: 1000
});

// Combine the SWIR1 bands from each image into a single image.
var swir1Bands = s2ImgSwir1.select('B11')
  .addBands(l8ImgSwir1.select('B6'))
  .rename(['s2_swir1', 'l8_swir1']);

// Sample the SWIR1 bands using the sample point feature collection.
var imgSamp = swir1Bands.sampleRegions({
  collection: sample,
  scale: 30
})
// Add a constant property to each feature to be used as an independent variable.
.map(function(feature) {
  return feature.set('constant', 1);
});

// Compute linear regression coefficients. numX is 2 because
// there are two independent variables: 'constant' and 's2_swir1'. numY is 1
// because there is a single dependent variable: 'l8_swir1'. Cast the resulting
// object to an ee.Dictionary for easy access to the properties.
var linearRegression = ee.Dictionary(imgSamp.reduceColumns({
  reducer: ee.Reducer.linearRegression({
    numX: 2,
    numY: 1
  }),
  selectors: ['constant', 's2_swir1', 'l8_swir1']
}));

// Convert the coefficients array to a list.
var coefList = ee.Array(linearRegression.get('coefficients')).toList();

// Extract the y-intercept and slope.
var yInt = ee.List(coefList.get(0)).get(0); // y-intercept
var slope = ee.List(coefList.get(1)).get(0); // slope

// Gather the SWIR1 values from the point sample into a list of lists.
var props = ee.List(['s2_swir1', 'l8_swir1']);
var regressionVarsList = ee.List(imgSamp.reduceColumns({
  reducer: ee.Reducer.toList().repeat(props.size()),
  selectors: props
}).get('list'));

// Convert regression x and y variable lists to an array - used later as input
// to ui.Chart.array.values for generating a scatter plot.
var x = ee.Array(ee.List(regressionVarsList.get(0)));
var y1 = ee.Array(ee.List(regressionVarsList.get(1)));

// Apply the line function defined by the slope and y-intercept of the
// regression to the x variable list to create an array that will represent
// the regression line in the scatter plot.
var y2 = ee.Array(ee.List(regressionVarsList.get(0)).map(function(x) {
  var y = ee.Number(x).multiply(slope).add(yInt);
  return y;
}));

// Concatenate the y variables (Landsat 8 SWIR1 and predicted y) into an array
// for input to ui.Chart.array.values for plotting a scatter plot.
var yArr = ee.Array.cat([y1, y2], 1);

// Make a scatter plot of the two SWIR1 bands for the point sample and include
// the least squares line of best fit through the data.
print(ui.Chart.array.values({
  array: yArr,
  axis: 0,
  xLabels: x})
  .setChartType('ScatterChart')
  .setOptions({
    legend: {position: 'none'},
    hAxis: {'title': 'Sentinel-2 SWIR1'},
    vAxis: {'title': 'Landsat 8 SWIR1'},
    series: {
      0: {
        pointSize: 0.2,
        dataOpacity: 0.5,
      },
      1: {
        pointSize: 0,
        lineWidth: 2,
      }
    }
  })
);


Figura 2. Grafico a dispersione e linea di regressione lineare dei minimi quadrati per un campione di pixel che rappresentano la riflettanza TOA di Sentinel-2 e Landsat 8 SWIR1.

ee.List

Le colonne di oggetti ee.List 2D possono essere inserite nei riduttori di regressione. Gli esempi seguenti forniscono prove semplici; la variabile indipendente è una copia della variabile dipendente che produce un'intercetta sull'asse y pari a 0 e una pendenza pari a 1.

linearFit()

Editor di codice (JavaScript)

// Define a list of lists, where columns represent variables. The first column
// is the independent variable and the second is the dependent variable.
var listsVarColumns = ee.List([
  [1, 1],
  [2, 2],
  [3, 3],
  [4, 4],
  [5, 5]
]);

// Compute the least squares estimate of a linear function. Note that an
// object is returned; cast it as an ee.Dictionary to make accessing the
// coefficients easier.
var linearFit = ee.Dictionary(listsVarColumns.reduce(ee.Reducer.linearFit()));

// Inspect the result.
print(linearFit);
print('y-intercept:', linearFit.get('offset'));
print('Slope:', linearFit.get('scale'));

Trasponi l'elenco se le variabili sono rappresentate da righe convertendole in un ee.Array, trasponendole e poi convertendole di nuovo in un ee.List.

Editor di codice (JavaScript)

// If variables in the list are arranged as rows, you'll need to transpose it.
// Define a list of lists where rows represent variables. The first row is the
// independent variable and the second is the dependent variable.
var listsVarRows = ee.List([
  [1, 2, 3, 4, 5],
  [1, 2, 3, 4, 5]
]);

// Cast the ee.List as an ee.Array, transpose it, and cast back to ee.List.
var listsVarColumns = ee.Array(listsVarRows).transpose().toList();

// Compute the least squares estimate of a linear function. Note that an
// object is returned; cast it as an ee.Dictionary to make accessing the
// coefficients easier.
var linearFit = ee.Dictionary(listsVarColumns.reduce(ee.Reducer.linearFit()));

// Inspect the result.
print(linearFit);
print('y-intercept:', linearFit.get('offset'));
print('Slope:', linearFit.get('scale'));

linearRegression()

L'applicazione di ee.Reducer.linearRegression() è simile all'esempio di linearFit() riportato sopra, tranne per il fatto che è inclusa una variabile indipendente costante.

Editor di codice (JavaScript)

// Define a list of lists where columns represent variables. The first column
// represents a constant term, the second an independent variable, and the third
// a dependent variable.
var listsVarColumns = ee.List([
  [1, 1, 1],
  [1, 2, 2],
  [1, 3, 3],
  [1, 4, 4],
  [1, 5, 5]
]);

// Compute ordinary least squares regression coefficients. numX is 2 because
// there is one constant term and an additional independent variable. numY is 1
// because there is only a single dependent variable. Cast the resulting
// object to an ee.Dictionary for easy access to the properties.
var linearRegression = ee.Dictionary(
  listsVarColumns.reduce(ee.Reducer.linearRegression({
    numX: 2,
    numY: 1
})));

// Convert the coefficients array to a list.
var coefList = ee.Array(linearRegression.get('coefficients')).toList();

// Extract the y-intercept and slope.
var b0 = ee.List(coefList.get(0)).get(0); // y-intercept
var b1 = ee.List(coefList.get(1)).get(0); // slope

// Extract the residuals.
var residuals = ee.Array(linearRegression.get('residuals')).toList().get(0);

// Inspect the results.
print('OLS estimates', linearRegression);
print('y-intercept:', b0);
print('Slope:', b1);
print('Residuals:', residuals);