采样点的选取
如果你采用监督学习的话，那就手动打标签 或者可以了解一下非监督学习
合成多季节多波段影像
首先，制作一个包含多波段的影像，每个波段作为随机森林分类器的一个feature输入，提升feature的丰富度以保证分类精度。
1、对landsat5用的云掩膜函数：
代码语言：javascript复制// cloud mask
var cloudMaskL457 = function(image) {
var qa = image.select('pixel_qa');
// If the cloud bit (5) is set and the cloud confidence (7) is high
// or the cloud shadow bit is set (3), then it's a bad pixel.
var cloud = qa.bitwiseAnd(1 << 5)
.and(qa.bitwiseAnd(1 << 7))
.or(qa.bitwiseAnd(1 << 3));
// Remove edge pixels that don't occur in all bands
var mask2 = image.mask().reduce(ee.Reducer.min());
return image.updateMask(cloud.not()).updateMask(mask2);
};2、构建正向指标和反指标
代码语言：javascript复制// indicators
// landsat 5
function NDWI(img){
var nir = img.select('B4');
var green = img.select('B2');
var ndwi = img.expression(
'(B2-B4)/(B2 B4)',
{
'B4':nir,
'B2':green
});
return ndwi;
}
function MNDWI(img){
var swir1 = img.select('B5');
var green = img.select('B2');
var mndwi = img.expression(
'(B2-B5)/(B2 B5)',
{
'B5':swir1,
'B2':green
});
return mndwi;
}
function EWI(img){
var swir1 = img.select('B5')
var nir = img.select('B4');
var green = img.select('B2');
var ewi = img.expression(
'(B2-B4-B5)/(B2 B4 B5)',
{
'B5':swir1,
'B4':nir,
'B2':green
});
return ewi;
}
function NDBI(img){
var swir1 = img.select('B5');
var nir = img.select('B4');
var ndbi = img.expression(
'(B5-B4)/(B5 B4)',
{
'B5':swir1,
'B4':nir
});
return ndbi;
}
function NDVI(img){
var nir = img.select('B4');
var red = img.select('B3');
var ndvi = img.expression(
'(B4-B3)/(B4 B3)',
{
'B4':nir,
'B3':red
});
return ndvi;
}3、选取我们的影像数据源，这里用的是landsat5，如果更换数据源的话云掩膜函数和指标计算也要一并调整。影像源代码如下：
代码语言：javascript复制// satellite data
var l8_col = ee.ImageCollection("LANDSAT/LT05/C01/T1_SR");4、开始合成季节影像。 这里我做的分类是1990年的，但是因为landsat的重返周期太长了，尺度稍微大点的话会碰上很多云，掩膜的话又都掩膜没了，所以我使用1989到1991三年的合成，其中冬是12月1日到3月1日、春是3月1日到6月1日、夏是6月1日到9月1日、秋是9月1日到12月1日，可以根据自己的需求调整。另外，为了避免云掩膜出现的孔洞，所以我用unmask函数把孔洞赋值为0以填补出现的孔洞。具体代码如下：
代码语言：javascript复制// winter
var img_winter = ee.List([])
for (var i=1989; i<1991; i ){
var winter = ee.Image(l8_col.filterBounds(roi)
.filterDate(i "-12-1", (i 1) "-3-1")
.filter(ee.Filter.lt('CLOUD_COVER', 100))
.map(cloudMaskL457)
.mean());
}
img_winter = img_winter.add(winter);
img_winter = ee.ImageCollection(img_winter).mean().unmask(0);
// spring
var img_spring = ee.List([])
for (var i=1989; i<1991; i ){
var spring = ee.Image(l8_col.filterBounds(roi)
.filterDate(i "-3-1", i "-6-1")
.filter(ee.Filter.lt('CLOUD_COVER', 100))
.map(cloudMaskL457)
.mean());
}
img_spring = img_spring.add(spring);
img_spring = ee.ImageCollection(img_spring).mean().unmask(0);
// summer
var img_summer = ee.List([])
for (var i=1989; i<1991; i ){
var summer = ee.Image(l8_col.filterBounds(roi)
.filterDate(i "-6-1", i "-9-1")
.filter(ee.Filter.lt('CLOUD_COVER', 100))
.map(cloudMaskL457)
.mean());
}
img_summer = img_summer.add(summer);
img_summer = ee.ImageCollection(img_summer).mean().unmask(0);
// fall
var img_fall = ee.List([])
for (var i=1989; i<1991; i ){
var fall = ee.Image(l8_col.filterBounds(roi)
.filterDate(i "-9-1", i "-12-1")
.filter(ee.Filter.lt('CLOUD_COVER', 100))
.map(cloudMaskL457)
.mean());
}
img_fall = img_fall.add(fall);
img_fall = ee.ImageCollection(img_fall).mean().unmask(0);5、为了确保分类成功率，我引入GEE中现成的几个产品。 这些产品可以根据需要灵活调整，这里就不多介绍了。依然地，使用unmask函数填补掩膜孔洞，代码如下：
代码语言：javascript复制// referring image data
var impervious = ee.Image("Tsinghua/FROM-GLC/GAIA/v10").clip(roi).unmask(-1);
var water = ee.Image('JRC/GSW1_2/GlobalSurfaceWater').clip(roi).unmask(-1);6、之前合成了landsat5季节影像，别忘记上面我们写下了各种指数 这里就要使用这些合成影像计算指数生成各指数的波段，并进行命名，代码如下：
代码语言：javascript复制var ndwi_wi = NDWI(img_winter).rename('ndwi_wi');
var mndwi_wi = MNDWI(img_winter).rename('mndwi_wi');
var ewi_wi = EWI(img_winter).rename('ewi_wi');
var ndbi_wi = NDBI(img_winter).rename('ndbi_wi');
var ndvi_wi = NDVI(img_spring).rename('ndvi_wi');
var ndwi_sp = NDWI(img_spring).rename('ndwi_sp');
var mndwi_sp = MNDWI(img_spring).rename('mndwi_sp');
var ewi_sp = EWI(img_spring).rename('ewi_sp');
var ndbi_sp = NDBI(img_spring).rename('ndbi_sp');
var ndvi_sp = NDVI(img_spring).rename('ndvi_sp');
var ndwi_su = NDWI(img_summer).rename('ndwi_su');
var mndwi_su = MNDWI(img_summer).rename('mndwi_su');
var ewi_su = EWI(img_summer).rename('ewi_su');
var ndbi_su = NDBI(img_summer).rename('ndbi_su');
var ndvi_su = NDVI(img_summer).rename('ndvi_su');
var ndwi_fa = NDWI(img_fall).rename('ndwi_fa');
var mndwi_fa = MNDWI(img_fall).rename('mndwi_fa');
var ewi_fa = EWI(img_fall).rename('ewi_fa');
var ndbi_fa = NDBI(img_fall).rename('ndbi_fa');
var ndvi_fa = NDVI(img_fall).rename('ndvi_fa');7、上面的指数是作为波段的形式。然后我们把引用的产品中的波段也提出并进行重命名，代码如下：
代码语言：javascript复制var imperchange = impervious.select('change_year_index').rename('imperchange');
var waterocc = water.select('occurrence').unmask(-1).rename('waterocc');
var waterchange = water.select('change_norm').unmask(-150).rename('waterchange');8、然后，我们从每幅季节合成影像中提取波段1-7，并且加入同一幅影像，代码如下：
代码语言：javascript复制// add bands into images
var image = img_winter.addBands([ndwi_wi, mndwi_wi, ewi_wi, ndbi_wi, ndvi_wi,
img_spring.select('B1').rename('B1_sp'), img_spring.select('B2').rename('B2_sp'), img_spring.select('B3').rename('B3_sp'), img_spring.select('B4').rename('B4_sp'), img_spring.select('B5').rename('B5_sp'),
ndwi_sp, mndwi_sp, ewi_sp, ndbi_sp, ndvi_sp,
img_summer.select('B1').rename('B1_su'), img_summer.select('B2').rename('B2_su'), img_summer.select('B3').rename('B3_su'), img_summer.select('B4').rename('B4_su'), img_summer.select('B5').rename('B5_su'),
ndwi_su, mndwi_su, ewi_su, ndbi_su, ndvi_su,
img_fall.select('B1').rename('B1_fa'), img_fall.select('B2').rename('B2_fa'), img_fall.select('B3').rename('B3_fa'), img_fall.select('B4').rename('B4_fa'), img_fall.select('B5').rename('B5_fa'),
ndwi_fa, mndwi_fa, ewi_fa, ndbi_fa, ndvi_fa,
imperchange, waterocc, waterchange]);分类器参数设置
上面我们合成了一个含有51波段的影像，我们现在就需要选取我们分类器需要输入的波段了，代码如下：
代码语言：javascript复制// select bands
var bands = ['B1', 'B2', 'B3', 'B4', 'B5', 'B6', 'B7',
'ndwi_wi', 'mndwi_wi', 'ewi_wi', 'ndbi_wi', 'ndvi_wi',
'B1_sp', 'B2_sp', 'B3_sp', 'B4_sp', 'B5_sp',
'ndwi_sp', 'mndwi_sp', 'ewi_sp', 'ndbi_sp', 'ndvi_sp',
'B1_su', 'B2_su', 'B3_su', 'B4_su', 'B5_su',
'ndwi_su', 'mndwi_su', 'ewi_su', 'ndbi_su', 'ndvi_su',
'B1_fa', 'B2_fa', 'B3_fa', 'B4_fa', 'B5_fa',
'ndwi_fa', 'mndwi_fa', 'ewi_fa', 'ndbi_fa', 'ndvi_fa',
'imperchange', 'waterocc', 'waterchange'];用merge函数把各分类合成一整个集合，然后它们的波段是value，然后把训练集和测试集按90%和10%分开。代码如下：
代码语言：javascript复制// set training data
var SP = Water.merge(Ponds).merge(Mud).merge(Others);
var training = image.select(bands).sampleRegions({
collection: SP,
properties:['value'],
scale: 30
});
// set varified data
var withRandom = training.randomColumn('random');
var split = 0.9;
var trainingP = withRandom.filter(ee.Filter.lt('random', split)); // 90% training data
var testingP = withRandom.filter(ee.Filter.gte('random', split)); // 10% verified data监督分类的实现与分类精度计算
得到一个分类后的影像
代码语言：javascript复制// classification
var class_img = image.select(bands).classify(classifier);光有分类结果不行，还需要知道分类精度。把测试集也分类一下，然后用GEE自带的函数求转移矩阵并且计算overall accuracy和kappa，这些会打印在console里，代码如下：
代码语言：javascript复制// testing
var test = testingP.classify(classifier);
// confusion matrix
var confusionMatrix = test.errorMatrix('value', 'classification')
print('confusionMatrix',confusionMatrix);//面板上显示混淆矩阵
print('overall accuracy', confusionMatrix.accuracy());//面板上显示总体精度
print('kappa accuracy', confusionMatrix.kappa());//面板上显示kappa值最后，当然要显示季节影像和分类后的影像，RGB波段组合还是色带什么的可以根据需要自行调整。如果选择手点的话，可以先随便点几个点，然后再根据显示的季节影像和分类影像再增加数据集，慢慢达到精度，代码如下：
代码语言：javascript复制// show images
var class_color = {
min: 0,
max: 3,
palette: ['ffffff', '0049ff', '00ffca', '97422c']
};
Map.addLayer(img_winter.clip(roi.geometry()), {bands: ["B4", "B3", "B2"], min:0, max:2000}, "winter");
Map.addLayer(img_spring.clip(roi.geometry()), {bands: ["B4", "B3", "B2"], min:0, max:2000}, "spring");
Map.addLayer(img_summer.clip(roi.geometry()), {bands: ["B4", "B3", "B2"], min:0, max:2000}, "summer");
Map.addLayer(img_fall.clip(roi.geometry()), {bands: ["B4", "B3", "B2"], min:0, max:2000}, "fall");
Map.addLayer(class_img.clip(roi.geometry()), class_color, 'classifed')
Map.centerObject(roi, 7)最后是保存在云盘里，下载
代码语言：javascript复制Export.image.toDrive({
image: class_img,//分类结果
description: 'test_season_1990',//文件名
folder: '111',
scale: 30,//分辨率
region: roi,//区域
crs:"EPSG:4326",
maxPixels:1e13//此处值设置大一些，防止溢出
});
