classification_regressionjs

/**** Start of imports. If edited, may not auto-convert in the playground. ****/
var forest = /* color: #589400 */ee.FeatureCollection(
        [ee.Feature(
            ee.Geometry.Point([9.239888569645757, 44.70488019367973]),
            {
              "class": 0,
              "system:index": "0"
            }),
        ee.Feature(
            ee.Geometry.Point([9.415669819645757, 44.65019546327906]),
            {
              "class": 0,
              "system:index": "1"
            }),
        ee.Feature(
            ee.Geometry.Point([9.585957905583257, 44.66973164857281]),
            {
              "class": 0,
              "system:index": "2"
            }),
        ee.Feature(
            ee.Geometry.Point([9.514546772770757, 44.57198489627976]),
            {
              "class": 0,
              "system:index": "3"
            }),
        ee.Feature(
            ee.Geometry.Point([9.382710835270757, 44.60719268481843]),
            {
              "class": 0,
              "system:index": "4"
            }),
        ee.Feature(
            ee.Geometry.Point([9.585957905583257, 44.60719268481843]),
            {
              "class": 0,
              "system:index": "5"
            }),
        ee.Feature(
            ee.Geometry.Point([9.657369038395757, 44.63847059243192]),
            {
              "class": 0,
              "system:index": "6"
            }),
        ee.Feature(
            ee.Geometry.Point([9.800191304020757, 44.61501374153445]),
            {
              "class": 0,
              "system:index": "7"
            })]),
    developed = /* color: #ff0000 */ee.FeatureCollection(
        [ee.Feature(
            ee.Geometry.Point([8.899312397770757, 44.46623357055405]),
            {
              "class": 1,
              "system:index": "0"
            }),
        ee.Feature(
            ee.Geometry.Point([9.086079975895757, 44.423094670984064]),
            {
              "class": 1,
              "system:index": "1"
            }),
        ee.Feature(
            ee.Geometry.Point([10.360494038395757, 44.845261089622966]),
            {
              "class": 1,
              "system:index": "2"
            }),
        ee.Feature(
            ee.Geometry.Point([8.937764546208257, 44.509340613304694]),
            {
              "class": 1,
              "system:index": "3"
            }),
        ee.Feature(
            ee.Geometry.Point([9.198689839177007, 45.21407215238483]),
            {
              "class": 1,
              "system:index": "4"
            }),
        ee.Feature(
            ee.Geometry.Point([9.297566792302007, 45.19084912934457]),
            {
              "class": 1,
              "system:index": "5"
            })]),
    water = /* color: #1a11ff */ee.FeatureCollection(
        [ee.Feature(
            ee.Geometry.Point([10.371480366520757, 44.36421756103338]),
            {
              "class": 2,
              "system:index": "0"
            }),
        ee.Feature(
            ee.Geometry.Point([8.866353413395757, 44.17147289167903]),
            {
              "class": 2,
              "system:index": "1"
            }),
        ee.Feature(
            ee.Geometry.Point([9.031148335270757, 44.202984579473394]),
            {
              "class": 2,
              "system:index": "2"
            }),
        ee.Feature(
            ee.Geometry.Point([9.195943257145757, 44.175412773917195]),
            {
              "class": 2,
              "system:index": "3"
            }),
        ee.Feature(
            ee.Geometry.Point([9.025655171208257, 44.15965166567795]),
            {
              "class": 2,
              "system:index": "4"
            }),
        ee.Feature(
            ee.Geometry.Point([8.849873921208257, 44.317073224769544]),
            {
              "class": 2,
              "system:index": "5"
            }),
        ee.Feature(
            ee.Geometry.Point([8.761983296208257, 44.21085986920881]),
            {
              "class": 2,
              "system:index": "6"
            }),
        ee.Feature(
            ee.Geometry.Point([8.948750874333257, 44.195108236844554]),
            {
              "class": 2,
              "system:index": "7"
            })]),
    herbaceous = /* color: #ffc82d */ee.FeatureCollection(
        [ee.Feature(
            ee.Geometry.Point([10.942769429020757, 44.79460735427151]),
            {
              "class": 3,
              "system:index": "0"
            }),
        ee.Feature(
            ee.Geometry.Point([10.871358296208257, 44.85694409346852]),
            {
              "class": 3,
              "system:index": "1"
            }),
        ee.Feature(
            ee.Geometry.Point([10.756001850895757, 44.868624727416375]),
            {
              "class": 3,
              "system:index": "2"
            }),
        ee.Feature(
            ee.Geometry.Point([10.569234272770757, 44.923103017066865]),
            {
              "class": 3,
              "system:index": "3"
            }),
        ee.Feature(
            ee.Geometry.Point([10.201192280583257, 44.88419521949802]),
            {
              "class": 3,
              "system:index": "4"
            }),
        ee.Feature(
            ee.Geometry.Point([9.986958882145757, 44.930881416692465]),
            {
              "class": 3,
              "system:index": "5"
            }),
        ee.Feature(
            ee.Geometry.Point([10.750508686833257, 44.84915568755982]),
            {
              "class": 3,
              "system:index": "6"
            })]);
/***** End of imports. If edited, may not auto-convert in the playground. *****/
//the processes of training data collection, 
//assifier selection, classifier training, and image classification

// Create an Earth Engine Point object over Milan. 
var pt = ee.Geometry.Point([9.453, 45.424]); 
// Filter the Landsat 8 collection and select the least cloudy image. 
var landsat = ee.ImageCollection('LANDSAT/LC08/C02/T1_L2') 
.filterBounds(pt) 
.filterDate('2019-01-01', '2020-01-01') 
.sort('CLOUD_COVER') 
.first(); 
// Center the map on that image. 
Map.centerObject(landsat, 8); 
// Add Landsat image to the map. 
var visParams = { 
bands: ['SR_B4', 'SR_B3', 'SR_B2'], 
min: 7000, 
max: 12000 
}; 
Map.addLayer(landsat, visParams, 'Landsat 8 image'
);

// Combine training feature collections. 
var trainingFeatures = ee.FeatureCollection([ 
forest, developed, water, herbaceous 
]).flatten();

// Define prediction bands. 
var predictionBands = [ 
'SR_B1', 'SR_B2', 'SR_B3', 'SR_B4', 'SR_B5', 'SR_B6', 
'SR_B7','ST_B10' 
]; 
// Sample training points. 
var classifierTraining = landsat.select(predictionBands) 
.sampleRegions({ 
collection: trainingFeatures, 
properties: ['class'], 
scale: 30 
});

//////////////// CART Classifier /////////////////// 
// Train a CART Classifier. 
var classifier = ee.Classifier.smileCart().train({ 
features: classifierTraining, 
classProperty: 'class', 
inputProperties: predictionBands }
);

// Classify the Landsat image. 
var classified = 
landsat.select(predictionBands).classify(classifier); 
// Define classification image visualization parameters. 
var classificationVis = { 
min: 0, 
max: 3, 
palette: ['589400', 'ff0000', '1a11ff', 'd0741e'] 
}; 
// Add the classified image to the map. 
Map.addLayer(classified, classificationVis, 'CART classified');

 /////////////// Random Forest Classifier 
///////////////////// 
// Train RF classifier. 
var RFclassifier = 
ee.Classifier.smileRandomForest(50).train({ 
features: classifierTraining, 
classProperty: 'class', 
inputProperties: predictionBands 
}); 
// Classify Landsat image. 
var RFclassified = 
landsat.select(predictionBands).classify( 
RFclassifier); 
// Add classified image to the map. 
Map.addLayer(RFclassified, classificationVis, 'RF classified');
//////////////////////////////////////////
///////////////////// Unsupervised Classification
//////////////////////////////////////////

// Make the training dataset. 
var training = landsat.sample({ 
region: landsat.geometry(), 
scale: 30, 
numPixels: 1000, 
tileScale: 8});

 // Instantiate the clusterer and train it. 
var clusterer = 
ee.Clusterer.wekaKMeans(4).train(training);

 // Cluster the input using the trained clusterer. 
var Kclassified = landsat.cluster(clusterer); 
// Display the clusters with random colors. 
Map.addLayer(Kclassified.randomVisualizer(), {}, 
'K-means classified -random colors');


////////////////////////////////////////////////////////////////////////////////////////
// Import the reference dataset. 
var data = ee.FeatureCollection( 
'projects/gee-book/assets/F2-2/milan_data'); 
// Define the prediction bands. 
var predictionBands = [ 
'SR_B1', 'SR_B2', 'SR_B3', 'SR_B4', 'SR_B5', 'SR_B6', 
'SR_B7','ST_B10','ndvi', 'ndwi' 
]; 
// Split the dataset into training and testing sets. 
var trainingTesting = data.randomColumn(); 
var trainingSet = trainingTesting 
.filter(ee.Filter.lessThan('random', 0.8)); 
var testingSet = trainingTesting 
.filter(ee.Filter.greaterThanOrEquals('random', 0.8));

// Train the Random Forest Classifier with the trainingSet. 
var RFclassifier = 
ee.Classifier.smileRandomForest(50).train({ 
features: trainingSet, 
classProperty: 'class', 
inputProperties: predictionBands 
});
// Now, to test the classification (verify model's accuracy), 
// we classify the testingSet and get a confusion matrix. 
var confusionMatrix = testingSet.classify(RFclassifier) 
.errorMatrix({ 
actual: 'class', 
predicted: 'classification' 
});

// Print the results. 
print('Confusion matrix:', confusionMatrix); 
print('Overall Accuracy:', confusionMatrix.accuracy()); 
print('Producers Accuracy:', 
confusionMatrix.producersAccuracy()); 
print('Consumers Accuracy:', 
confusionMatrix.consumersAccuracy()); 
print('Kappa:', confusionMatrix.kappa());

 // Hyperparameter tuning. 
var numTrees = ee.List.sequence(5, 100, 5); 
var accuracies = numTrees.map(function(t) { 
var classifier = ee.Classifier.smileRandomForest(t) 
.train({ 
features: trainingSet, 
classProperty: 'class', 
inputProperties: predictionBands 
}); 
return testingSet 
.classify(classifier) 
.errorMatrix('class', 'classification') 
.accuracy(); 
}); 
print(ui.Chart.array.values({ 
array: ee.Array(accuracies), 
axis: 0, 
xLabels: numTrees 
}).setOptions({ 
hAxis: { 
title: 'Number of trees' 
}, 
vAxis: { 
title: 'Accuracy' 
}, 
title: 'Accuracy per number of trees' 
})); 


//////////////////////////////ch8
///////Interpreting an Image: Regression
//////////////////////////////

// Define a Turin polygon. 
var Turin = ee.Geometry.Polygon( 
[ 
[ 
[7.455553918110218, 45.258245019259036], 
[7.455553918110218, 44.71237367431335], 
[8.573412804828967, 44.71237367431335], 
[8.573412804828967, 45.258245019259036] 
] 
], null, false); 
// Center on Turin 
Map.centerObject(Turin, 9); 
var mod44b = ee.ImageCollection('MODIS/006/MOD44B');
///// 
// Start Linear Fit 
///// 
// Put together the dependent variable by filtering the 
// ImageCollection to just the 2020 image near Turin and 
// selecting the percent tree cover band. 
var percentTree2020 = mod44b 
.filterDate('2020-01-01', '2021-01-01') 
.first() 
.clip(Turin) 
.select('Percent_Tree_Cover'); 
// You can print information to the console for inspection. 
print('2020 Image', percentTree2020); 
Map.addLayer(percentTree2020, { 
max: 100 
}, 'Percent Tree Cover');
var landsat8_raw = 
ee.ImageCollection('LANDSAT/LC08/C02/T1_RT');
// Put together the independent variable. 
var landsat8filtered = landsat8_raw 
.filterBounds(Turin.centroid({ 
'maxError': 1 
})) 
.filterDate('2020-04-01', '2020-4-30') 
.first(); 
print('Landsat8 filtered', landsat8filtered); 
// Display the L8 image. 
var visParams = { 
bands: ['B4', 'B3', 'B2'], 
max: 16000 
}; 
Map.addLayer(landsat8filtered, visParams, 'Landsat 8 Image');
// Calculate NDVI which will be the independent variable. 
var ndvi = landsat8filtered.normalizedDifference(['B5', 
'B4']); 
// Create the training image. 
var trainingImage = ndvi.addBands(percentTree2020); 
print('training image for linear fit', trainingImage);

// Independent variable first, dependent variable second. 
// You need to include the scale variable. 
var linearFit = trainingImage.reduceRegion({ 
reducer: ee.Reducer.linearFit(), 
geometry: Turin, 
scale: 30, 
bestEffort: true 
}); 
// Inspect the results. 
print('OLS estimates:', linearFit); 
print('y-intercept:', linearFit.get('offset')); 
print('Slope:', linearFit.get('scale'));
// Create a prediction based on the linearFit model. 
var predictedTree = ndvi.expression( 
'intercept + slope * ndvi', { 
'ndvi': ndvi.select('nd'), 
'intercept': ee.Number(linearFit.get('offset')), 
'slope': ee.Number(linearFit.get('scale')) 
}); 
print('predictedTree', predictedTree); 
// Display the results. 
Map.addLayer(predictedTree, { 
max: 100 
}, 'Predicted Percent Tree Cover');

////// 
// Start Linear Regression 
////// 
// Assemble the independent variables. 
var predictionBands = ['B1', 'B2', 'B3', 'B4', 'B5', 'B6', 
'B7','B10', 'B11' ];

// Create the training image stack for linear regression. 
var trainingImageLR = ee.Image(1) 
.addBands(landsat8filtered.select(predictionBands)) 
.addBands(percentTree2020); 
print('Linear Regression training image:', trainingImageLR);

 // Compute ordinary least squares regression coefficients using 
// the linearRegression reducer. 
var linearRegression = trainingImageLR.reduceRegion({ 
reducer: ee.Reducer.linearRegression({ 
numX: 10, 
numY: 1 
}), 
geometry: Turin, 
scale: 30, 
bestEffort: true 
});

// Inspect the results. 
print('Linear regression results:', linearRegression);
// Extract the coefficients as a list. 
var coefficients = 
ee.Array(linearRegression.get('coefficients')) 
.project([0]) 
.toList(); 
print('Coefficients', coefficients);

// Create the predicted tree cover based on linear regression. 
var predictedTreeLR = ee.Image(1) 
.addBands(landsat8filtered.select(predictionBands)) 
.multiply(ee.Image.constant(coefficients)) 
.reduce(ee.Reducer.sum()) 
.rename('predictedTreeLR') 
.clip(landsat8filtered.geometry()); 
Map.addLayer(predictedTreeLR, { 
min: 0, 
max: 100 
}, 'LR prediction');

///// 
// Start Non-linear Regression 
///// 
// Create the training data stack. 
var trainingImageCART = 
ee.Image(landsat8filtered.select(predictionBands)) 
.addBands(percentTree2020);

// Sample the training data stack. 
var trainingData = trainingImageCART.sample({ 
region: Turin, 
scale: 30, 
numPixels: 1500, 
seed: 5 
}); 
// Examine the CART training data. 
print('CART training data',trainingData);

// Run the CART regression. 
var cartRegression = ee.Classifier.smileCart() 
.setOutputMode('REGRESSION') 
.train({ 
features: trainingData, 
classProperty: 'Percent_Tree_Cover', 
inputProperties: predictionBands });

// Run the CART regression. 
var cartRegression = ee.Classifier.smileCart() 
.setOutputMode('REGRESSION') 
.train({ 
features: trainingData, 
classProperty: 'Percent_Tree_Cover', 
inputProperties: predictionBands });

 // Create a prediction of tree cover based on the CART regression. 
var cartRegressionImage = 
landsat8filtered.select(predictionBands) 
.classify(cartRegression, 'cartRegression'); 
Map.addLayer(cartRegressionImage, { 
min: 0, 
max: 100 
}, 'CART regression');


///// 
// Calculating RMSE to assess model performance 
///// 
// Concatenate percent tree cover image and all predictions. 
var concat = ee.Image.cat(percentTree2020, 
predictedTree, 
predictedTreeLR, 
cartRegressionImage) 
.rename([ 
'TCpercent', 
'LFprediction', 
'LRprediction', 
'CARTprediction' 
]); 
// Sample pixels 
var sample = concat.sample({ 
region: Turin, 
scale: 30, 
numPixels: 500, 
seed: 5 
}); 
print('First feature in sample', sample.first());

// First step: This function computes the squared difference between 
// the predicted percent tree cover and the known percent tree cover 
var calculateDiff = function(feature) { 
var TCpercent = ee.Number(feature.get('TCpercent')); 
var diffLFsq = ee.Number(feature.get('LFprediction')) 
.subtract(TCpercent).pow(2); 
var diffLRsq = ee.Number(feature.get('LRprediction')) 
.subtract(TCpercent).pow(2); 
var diffCARTsq = 
ee.Number(feature.get('CARTprediction')) 
.subtract(TCpercent).pow(2); 
// Return the feature with the squared difference set to a 'diff' property. 
return feature.set({ 
'diffLFsq': diffLFsq, 
'diffLRsq': diffLRsq, 
'diffCARTsq': diffCARTsq 
});};
 // Second step: Calculate RMSE for population of pairs. 
var rmse = ee.List([ee.Number( 
// Map the difference function over the collection. 
sample.map(calculateDiff) 
// Reduce to get the mean squared difference. 
.reduceColumns({ 
reducer: ee.Reducer.mean().repeat(3), 
selectors: ['diffLFsq', 'diffLRsq', 
'diffCARTsq' 
] 
}).get('mean') 
// Flatten the list of lists. 
)]).flatten().map(function(i) { 
// Take the square root of the mean square differences. 
return ee.Number(i).sqrt(); 
}); 
// Print the result. 
print('RMSE', rmse);