في المقال التعليمي السابق (مقدمة)، رأينا كيف تسجّل عمليات التضمين المستندة إلى صور الأقمار الصناعية المسارات السنوية لبيانات الأقمار الصناعية ومتغيرات المناخ. ويجعل ذلك مجموعة البيانات قابلة للاستخدام بسهولة في تحديد مواقع المحاصيل بدون الحاجة إلى تصميم نموذج لعلم الأحياء الظاهري للمحاصيل. يُعدّ ربط أنواع المحاصيل مهمة صعبة تتطلّب عادةً وضع نماذج لظواهر نمو المحاصيل وجمع عيّنات من الحقول لجميع المحاصيل المزروعة في المنطقة.
في هذا البرنامج التعليمي، سنتّبع نهج التصنيف غير الخاضع للإشراف في رسم خرائط المحاصيل، ما يتيح لنا تنفيذ هذه المهمة المعقّدة بدون الاعتماد على تصنيفات الحقول. تستفيد هذه الطريقة من المعرفة المحلية بالمنطقة بالإضافة إلى إحصاءات المحاصيل المجمّعة، والتي تتوفّر بسهولة في العديد من مناطق العالم.
اختيار منطقة
في هذا البرنامج التعليمي، سننشئ خريطة لأنواع المحاصيل في مقاطعة سيرو غوردو، أيوا. تقع هذه المقاطعة في حزام الذرة في الولايات المتحدة الذي يضم محصولَين أساسيَّين، وهما الذرة وفول الصويا. هذه المعرفة المحلية مهمة وستساعدنا في تحديد المَعلمات الرئيسية لنموذجنا.
لنبدأ بالحصول على حدود المقاطعة المحدّدة.
// 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
بعد ذلك، نحمّل مجموعة بيانات "تضمين صور الأقمار الصناعية"، ونفلتر الصور حسب السنة المحدّدة، وننشئ صورة مجمّعة.
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) للحصول على قناع للمحاصيل.
// 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);
علينا أخذ صورة "تضمين القمر الصناعي" والحصول على عيّنات عشوائية لتدريب نموذج تجميع. بما أنّ المنطقة التي تهمّنا تحتوي على العديد من وحدات البكسل المخفية، قد يؤدي أخذ عيّنات عشوائية بسيطة إلى الحصول على عيّنات ذات قيم فارغة. لضمان إمكانية استخراج العدد المطلوب من العيّنات غير الفارغة، نستخدم أسلوب أخذ العيّنات الطبقي للحصول على العدد المطلوب من العيّنات في المناطق غير المحجوبة.
// 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);
الشكل: عيّنات عشوائية مستخرَجة للتجميع
تنفيذ التجميع غير الخاضع للإشراف
يمكننا الآن تدريب أداة تجميع وتقسيم متجهات التضمين 64D إلى عدد محدّد من المجموعات المميّزة. استنادًا إلى معرفتنا المحلية، هناك نوعان رئيسيان من المحاصيل يمثّلان غالبية المنطقة، بالإضافة إلى عدة محاصيل أخرى تغطي الجزء المتبقي. يمكننا إجراء تجميع عنقودي غير موجّه على "تضمين صور الأقمار الصناعية" للحصول على مجموعات من وحدات البكسل التي تتضمّن مسارات وأنماط زمنية متشابهة. سيتم تجميع وحدات البكسل التي تتضمّن خصائص طيفية ومكانية متشابهة بالإضافة إلى علم الظواهر الجوية المتشابهة في المجموعة نفسها.
تسمح لنا السمة 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);
نختار المجموعتين العنقوديتين اللتين تتضمّنان أكبر مساحة.
var selectedFc = clusterAreaFc.sort('area', false).limit(2);
print('Top 2 Clusters by Area', selectedFc);
لكننا ما زلنا لا نعرف أي مجموعة عنقودية تمثّل أي محصول. إذا كان لديك بعض عيّنات الحقول من الذرة أو فول الصويا، يمكنك تراكبها على المجموعات لتحديد تصنيفاتها. في حال عدم توفّر عيّنات من الحقول، يمكننا الاستفادة من إحصاءات المحاصيل المجمّعة. في العديد من مناطق العالم، يتم جمع إحصاءات المحاصيل الإجمالية ونشرها بانتظام. في الولايات المتحدة، تقدّم "الهيئة الوطنية للإحصاءات الزراعية" (NASS) إحصاءات تفصيلية عن المحاصيل لكل مقاطعة ولكل محصول رئيسي. في عام 2022، بلغت مساحة زراعة الذرة في مقاطعة سيرو غوردو بولاية آيوا 161,500 فدان، وبلغت مساحة زراعة فول الصويا 110,500 فدان.
باستخدام هذه المعلومات، نعرف الآن أنّ من بين المجموعتَين الرئيسيتَين، ستكون المجموعة التي تشغل أكبر مساحة هي على الأرجح حقول الذرة، بينما ستكون المجموعة الأخرى حقول فول الصويا. لنحدّد هذه التصنيفات ونقارن المساحات المحسوبة بالإحصاءات المنشورة.
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)');
للمساعدة في تفسير النتائج، يمكننا أيضًا استخدام عناصر واجهة المستخدم لإنشاء وسيلة إيضاح وإضافتها إلى الخريطة.
// 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'));
الشكل: خريطة المحاصيل التي تم رصدها والتي تتضمّن محاصيل الذرة وفول الصويا
التحقّق من صحة النتائج
تمكّنا من الحصول على خريطة لأنواع المحاصيل باستخدام مجموعة بيانات Satellite Embedding بدون أي تصنيفات للحقول، وذلك باستخدام الإحصاءات المجمّعة والمعرفة المحلية بالمنطقة فقط. لنقارن نتائجنا بخريطة أنواع المحاصيل الرسمية من "طبقات بيانات الأراضي الزراعية" (CDL) التابعة لهيئة NASS في وزارة الزراعة الأمريكية.
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)');
على الرغم من وجود اختلافات بين نتائجنا والخريطة الرسمية، ستلاحظ أنّنا تمكّنا من الحصول على نتائج جيدة جدًا بأقل جهد ممكن. من خلال تطبيق خطوات ما بعد المعالجة على النتائج، يمكننا إزالة بعض التشويش وملء الفجوات في الناتج.
جرِّب النص البرمجي الكامل لهذا الدليل التعليمي في "أداة تعديل الرموز البرمجية" في Earth Engine.