طبقه بندی بدون نظارت با مجموعه داده جاسازی ماهواره ای

در GitHub ویرایش کنید
گزارش مشکل
تاریخچه صفحه
نویسنده(ها): اندیشه های فضایی
این آموزش بخشی از مجموعه ای از آموزش های مربوط به مجموعه داده های جاسازی ماهواره است، همچنین به مقدمه ، طبقه بندی نظارت شده ، رگرسیون و جستجوی مشابه مراجعه کنید.

در آموزش قبلی (مقدمه)، دیدیم که چگونه Satellite Embeddings مسیرهای سالانه مشاهدات ماهواره ای و متغیرهای آب و هوایی را ثبت می کند. این باعث می شود مجموعه داده به راحتی برای نقشه برداری محصولات بدون نیاز به مدل سازی فنولوژی محصول قابل استفاده باشد. نقشه برداری از نوع محصول یک کار چالش برانگیز است که به طور معمول نیازمند مدل سازی فنولوژی محصول و جمع آوری نمونه های مزرعه برای همه محصولات کشت شده در منطقه است.

در این آموزش، ما یک رویکرد طبقه‌بندی بدون نظارت را برای نقشه‌برداری محصول خواهیم داشت که ما را قادر می‌سازد این کار پیچیده را بدون تکیه بر برچسب‌های زمینه انجام دهیم. این روش از دانش محلی منطقه به همراه آمار کل محصول استفاده می کند که به راحتی برای بسیاری از نقاط جهان در دسترس است.

یک منطقه را انتخاب کنید

برای این آموزش، ما یک نقشه نوع برش برای Cerro Gordo County، آیووا ایجاد خواهیم کرد. این شهرستان در کمربند ذرت ایالات متحده است که دارای دو محصول اصلی است: ذرت و سویا. این دانش محلی مهم است و به ما کمک می کند تا در مورد پارامترهای کلیدی مدل خود تصمیم گیری کنیم.

بیایید با به دست آوردن مرز برای شهرستان انتخاب شده شروع کنیم.

// Select the region
// Cerro Gordo County, Iowa
var counties = ee.FeatureCollection('TIGER/2018/Counties');

var selected = counties
  .filter(ee.Filter.eq('GEOID', '19033'));
var geometry = selected.geometry();
Map.centerObject(geometry, 12);
Map.addLayer(geometry, {color: 'red'}, 'Selected Region', false);


شکل: منطقه انتخاب شده

مجموعه داده های Satellite Embedding را آماده کنید

در مرحله بعد، مجموعه داده های Satellite Embedding را بارگذاری می کنیم، تصاویر را برای سال انتخابی فیلتر می کنیم و یک موزاییک ایجاد می کنیم.

var embeddings = ee.ImageCollection('GOOGLE/SATELLITE_EMBEDDING/V1/ANNUAL');

var year = 2022;
var startDate = ee.Date.fromYMD(year, 1, 1);
var endDate = startDate.advance(1, 'year');

var filteredembeddings = embeddings
  .filter(ee.Filter.date(startDate, endDate))
  .filter(ee.Filter.bounds(geometry));

var embeddingsImage = filteredembeddings.mosaic();

یک ماسک زمین زراعی ایجاد کنید

برای مدل سازی خود، باید مناطق غیر زراعی را حذف کنیم. مجموعه داده های جهانی و منطقه ای زیادی وجود دارد که می توان از آنها برای ایجاد یک ماسک برش استفاده کرد. ESA WorldCover یا GFSAD Global Cropland Extent Product انتخاب های خوبی برای مجموعه داده های جهانی زمین های زراعی هستند. محصول جدیدتر ESA WorldCereal Active Cropland است که نقشه‌برداری فصلی از زمین‌های زراعی فعال دارد. از آنجایی که منطقه ما در ایالات متحده است، می‌توانیم از مجموعه داده‌های منطقه‌ای دقیق‌تر USDA NASS Cropland Data Layers (CDL) برای بدست آوردن ماسک برش استفاده کنیم.

// Use Cropland Data Layers (CDL) to obtain cultivated cropland
var cdl = ee.ImageCollection('USDA/NASS/CDL')
  .filter(ee.Filter.date(startDate, endDate))
  .first();
var cropLandcover = cdl.select('cropland');
var croplandMask = cdl.select('cultivated').eq(2).rename('cropmask');

// Visualize the crop mask
var croplandMaskVis = {min: 0, max: 1, palette: ['white', 'green']};
Map.addLayer(croplandMask.clip(geometry), croplandMaskVis, 'Crop Mask');


شکل: منطقه انتخاب شده با ماسک زمین زراعی

نمونه های آموزشی را استخراج کنید

ماسک زراعی را روی موزاییک تعبیه شده اعمال می کنیم. ما اکنون با تمام پیکسل هایی که نشان دهنده زمین های زراعی کشت شده در شهرستان هستند باقی مانده ایم.

// Mask all non-cropland pixels
var clusterImage = embeddingsImage.updateMask(croplandMask);

برای آموزش یک مدل خوشه‌بندی باید تصویر Satellite Embedding را بگیریم و نمونه‌های تصادفی به دست آوریم. از آنجایی که منطقه مورد نظر ما حاوی بسیاری از پیکسل های ماسک شده است، یک نمونه گیری تصادفی ساده ممکن است به نمونه هایی با مقادیر تهی منجر شود. برای اطمینان از اینکه بتوانیم تعداد مورد نظر نمونه غیر پوچ را استخراج کنیم، از نمونه برداری طبقه بندی شده برای به دست آوردن تعداد مورد نظر نمونه در مناطق بدون نقاب استفاده می کنیم.

// Stratified random sampling
var training = clusterImage.addBands(croplandMask).stratifiedSample({
  numPoints: 1000,
  classBand: 'cropmask',
  region: geometry,
  scale: 10,
  tileScale: 16,
  seed: 100,
  dropNulls: true,
  geometries: true
});

صادرات نمونه به دارایی (اختیاری)

استخراج نمونه ها یک عملیات محاسباتی پرهزینه است و این روش خوبی است که نمونه های آموزشی استخراج شده را به عنوان دارایی صادر کنید و از دارایی های صادر شده در مراحل بعدی استفاده کنید. این به غلبه بر زمان محاسبات یا بیش از خطاهای حافظه کاربر هنگام کار با مناطق بزرگ کمک می کند.

کار صادرات را شروع کنید و منتظر بمانید تا قبل از ادامه کار به پایان برسد.

// Replace this with your asset folder
// The folder must exist before exporting
var exportFolder = 'projects/spatialthoughts/assets/satellite_embedding/';

var samplesExportFc = 'cluster_training_samples';
var samplesExportFcPath = exportFolder + samplesExportFc;

Export.table.toAsset({
  collection: training,
  description: 'Cluster_Training_Samples',
  assetId: samplesExportFcPath
});

پس از اتمام کار صادرات، می‌توانیم نمونه‌های استخراج‌شده را به عنوان مجموعه ویژگی‌ها در کد خود بازخوانی کنیم.

// Use the exported asset
var training = ee.FeatureCollection(samplesExportFcPath);

نمونه ها را تجسم کنید

چه نمونه خود را به صورت تعاملی اجرا کنید و چه به یک مجموعه ویژگی صادر کرده باشید، اکنون یک متغیر آموزشی با نقاط نمونه خود خواهید داشت. بیایید اولین نمونه را برای بررسی و اضافه کردن نقاط آموزشی خود به Map چاپ کنیم.

print('Extracted sample', training.first());
Map.addLayer(training, {color: 'blue'}, 'Extracted Samples', false);


شکل: نمونه های تصادفی استخراج شده برای خوشه بندی

انجام خوشه بندی بدون نظارت

اکنون می‌توانیم یک خوشه‌ساز را آموزش دهیم و بردارهای تعبیه‌شده 64 بعدی را در تعدادی از خوشه‌های مجزا گروه‌بندی کنیم. بر اساس دانش محلی ما، دو نوع اصلی از محصولات زراعی وجود دارد که اکثریت منطقه را با چندین محصول دیگر که بخش باقی مانده را پوشش می دهند، نشان می دهد. ما می‌توانیم خوشه‌بندی بدون نظارت را در Satellite Embedding انجام دهیم تا خوشه‌هایی از پیکسل‌ها را به‌دست آوریم که مسیرها و الگوهای زمانی مشابهی دارند. پیکسل هایی با ویژگی های طیفی و فضایی مشابه همراه با فنولوژی مشابه در یک خوشه گروه بندی می شوند.

ee.Clusterer.wekaCascadeKMeans() به ما اجازه می دهد حداقل و حداکثر تعداد خوشه ها را مشخص کنیم و تعداد بهینه خوشه ها را بر اساس داده های آموزشی پیدا کنیم. در اینجا دانش محلی ما برای تعیین حداقل و حداکثر تعداد خوشه ها مفید است. از آنجایی که انتظار داریم چند نوع متمایز از خوشه ها - ذرت، سویا، و چندین محصول دیگر - وجود داشته باشد، می توانیم از 4 به عنوان حداقل تعداد خوشه و 5 به عنوان حداکثر تعداد خوشه استفاده کنیم. ممکن است مجبور شوید این اعداد را آزمایش کنید تا بفهمید چه چیزی برای منطقه شما بهتر است.

// Cluster the Satellite Embedding Image
var minClusters = 4;
var maxClusters = 5;

var clusterer = ee.Clusterer.wekaCascadeKMeans({
  minClusters: minClusters, maxClusters: maxClusters}).train({
  features: training,
  inputProperties: clusterImage.bandNames()
});

var clustered = clusterImage.cluster(clusterer);
Map.addLayer(clustered.randomVisualizer().clip(geometry), {}, 'Clusters');


شکل: خوشه های به دست آمده از طبقه بندی بدون نظارت

برچسب ها را به خوشه ها اختصاص دهید

پس از بازرسی بصری، خوشه‌های به‌دست‌آمده در مراحل قبل با مرزهای مزرعه که در تصویر با وضوح بالا دیده می‌شوند، مطابقت دارند. ما از دانش محلی می دانیم که دو خوشه بزرگ ذرت و سویا هستند. بیایید مساحت هر خوشه را در تصویر خود محاسبه کنیم. 

// Calculate Cluster Areas
// 1 Acre = 4046.86 Sq. Meters
var areaImage = ee.Image.pixelArea().divide(4046.86).addBands(clustered);

var areas = areaImage.reduceRegion({
      reducer: ee.Reducer.sum().group({
      groupField: 1,
      groupName: 'cluster',
    }),
    geometry: geometry,
    scale: 10,
    maxPixels: 1e10
    });

var clusterAreas = ee.List(areas.get('groups'));

// Process results to extract the areas and create a FeatureCollection

var clusterAreas = clusterAreas.map(function(item) {
  var areaDict = ee.Dictionary(item);
  var clusterNumber = areaDict.getNumber('cluster').format();
  var area = areaDict.getNumber('sum');
  return ee.Feature(null, {cluster: clusterNumber, area: area});
});

var clusterAreaFc = ee.FeatureCollection(clusterAreas);
print('Cluster Areas', clusterAreaFc);

ما 2 خوشه را با بالاترین سطح انتخاب می کنیم. 

var selectedFc = clusterAreaFc.sort('area', false).limit(2);
print('Top 2 Clusters by Area', selectedFc);

اما هنوز نمی دانیم کدام خوشه چه محصولی است. اگر چند نمونه مزرعه ذرت یا سویا داشتید، می‌توانید آن‌ها را روی خوشه‌ها بپوشانید تا برچسب مربوط به آن‌ها را بفهمید. در غیاب نمونه‌های مزرعه، می‌توانیم از آمار کل محصول استفاده کنیم. در بسیاری از نقاط جهان، آمار کل محصولات به طور مرتب جمع آوری و منتشر می شود. برای ایالات متحده، خدمات ملی آمار کشاورزی (NASS) آمار دقیق محصول را برای هر شهرستان و هر محصول عمده دارد. برای سال 2022، شهرستان سرو گوردو، آیووا دارای مساحت زیر کشت ذرت 161500 هکتار و سطح زیر کشت سویا 110500 هکتار بود.

با استفاده از این اطلاعات، اکنون می دانیم که در میان 2 خوشه برتر، یکی با بیشترین سطح به احتمال زیاد ذرت و دیگری سویا خواهد بود. بیایید این برچسب ها را اختصاص دهیم و مناطق محاسبه شده را با آمار منتشر شده مقایسه کنیم. 

var cornFeature = selectedFc.sort('area', false).first();
var soybeanFeature = selectedFc.sort('area').first();
var cornCluster = cornFeature.get('cluster');
var soybeanCluster = soybeanFeature.get('cluster');

print('Corn Area (Detected)', cornFeature.getNumber('area').round());
print('Corn Area (From Crop Statistics)', 163500);

print('Soybean Area (Detected)', soybeanFeature.getNumber('area').round());
print('Soybean Area (From Crop Statistics)', 110500);

یک نقشه برش ایجاد کنید

اکنون برچسب‌های هر خوشه را می‌شناسیم و می‌توانیم پیکسل‌ها را برای هر نوع برش استخراج کنیم و آنها را برای ایجاد نقشه نهایی برش ادغام کنیم. 

// Select the clusters to create the crop map
var corn = clustered.eq(ee.Number.parse(cornCluster));
var soybean = clustered.eq(ee.Number.parse(soybeanCluster));

var merged = corn.add(soybean.multiply(2));
var cropVis = {min: 0, max: 2, palette: ['#bdbdbd', '#ffd400', '#267300']};
Map.addLayer(merged.clip(geometry), cropVis, 'Crop Map (Detected)');

برای کمک به تفسیر نتایج، می‌توانیم از عناصر UI برای ایجاد و اضافه کردن یک افسانه به نقشه استفاده کنیم. 

// Add a Legend
var legend = ui.Panel({
  layout: ui.Panel.Layout.Flow('horizontal'),
  style: {position: 'bottom-center', padding: '8px 15px'}});

var addItem = function(color, name) {
  var colorBox = ui.Label({
    style: {color: '#ffffff',
      backgroundColor: color,
      padding: '10px',
      margin: '0 4px 4px 0',
    }
  });
  var description = ui.Label({
    value: name,
    style: {
      margin: '0px 10px 0px 2px',
    }
  });
  return ui.Panel({
    widgets: [colorBox, description],
    layout: ui.Panel.Layout.Flow('horizontal')
  });
};

var title = ui.Label({
  value: 'Legend',
  style: {
    fontWeight: 'bold',
    fontSize: '16px',
    margin: '0px 10px 0px 4px'
  }
});

legend.add(title);
legend.add(addItem('#ffd400', 'Corn'));
legend.add(addItem('#267300', 'Soybean'));
legend.add(addItem('#bdbdbd', 'Other Crops'));


شکل: نقشه محصول شناسایی شده با محصولات ذرت و سویا

نتایج را تایید کنید

ما توانستیم یک نقشه نوع محصول را با مجموعه داده های جاسازی ماهواره ای بدون هیچ گونه برچسب میدانی با استفاده از آمار کل و دانش محلی منطقه به دست آوریم. بیایید نتایج خود را با نقشه رسمی نوع محصول از USDA NASS Cropland Data Layers (CDL) مقایسه کنیم. 

var cdl = ee.ImageCollection('USDA/NASS/CDL')
  .filter(ee.Filter.date(startDate, endDate))
  .first();
var cropLandcover = cdl.select('cropland');
var cropMap = cropLandcover.updateMask(croplandMask).rename('crops');

// Original data has unique values for each crop ranging from 0 to 254
var cropClasses = ee.List.sequence(0, 254);
// We remap all values as following
// Crop     | Source Value | Target Value
// Corn     | 1            | 1
// Soybean  | 5            | 2
// All other| 0-255        | 0
var targetClasses = ee.List.repeat(0, 255).set(1, 1).set(5, 2);
var cropMapReclass = cropMap.remap(cropClasses, targetClasses).rename('crops');

var cropVis = {min: 0, max: 2, palette: ['#bdbdbd', '#ffd400', '#267300']};
Map.addLayer(cropMapReclass.clip(geometry), cropVis, 'Crop Landcover (CDL)');


شکل: (سمت چپ) نقشه برش از جاسازی های ماهواره ای (راست) نقشه برش از CDL

در حالی که بین نتایج ما و نقشه رسمی اختلاف وجود دارد، متوجه خواهید شد که ما توانستیم با کمترین تلاش به نتایج بسیار خوبی دست یابیم. با اعمال مراحل پس از پردازش برای نتایج، می‌توانیم مقداری نویز را حذف کنیم و شکاف‌های خروجی را پر کنیم.

اسکریپت کامل این آموزش را در ویرایشگر کد موتور زمین امتحان کنید .