For nighttime light images

//Load china adm3 layer
var china = ee.FeatureCollection("projects/ee-chenyilin2025/assets/county");
var geometry = china.geometry();
Map.addLayer(china, {color: 'gray'}, 'Admin3 Boundaries');

//Load NTL image collection
var npp_viirs_ntl = ee.ImageCollection("projects/sat-io/open-datasets/npp-viirs-ntl");
var band = 'b1';

// Prepare MODIS Landcover collection
var modisLandcover = ee.ImageCollection('MODIS/061/MCD12Q1');
// Extract Landcover Classification for the year
// Using the IGBP classification
var classification = 'LC_Type1';

// Filter ntl image collection for the selected years
var startYear = 2001;
var endYear = 2020;
var startDate = ee.Date.fromYMD(startYear, 1, 1);
var endDate = ee.Date.fromYMD(endYear + 1, 1, 1);

var filtered = npp_viirs_ntl
  .filter(ee.Filter.date(startDate, endDate))
  .filter(ee.Filter.bounds(geometry))
  .select(band);
print('Filtered NTL collection', filtered);

var ntlImage = filtered.first();
var projection = ntlImage.projection();
var resolution = projection.nominalScale();

// Filter image collection for the selected years 
// and select the band having class values
var landcoverFiltered = modisLandcover
  .filter(ee.Filter.date(startDate, endDate))
  .select(classification);
  
print('Filtered MODIS Collection from ' + startYear + ' to ' + endYear,
  landcoverFiltered);

// Write a function to extract urban and agriculture classes from
// a MODIS landcover image
var extractClasses = function(image) {
  var urban = image.eq(13)
    .rename('urban');
  var agriculture = image.eq(10)
    .or(image.eq(12))
    .or(image.eq(14))
    .rename('agriculture')
 
  // Create a 2-band image with the extracted classes
  var classImage = urban.addBands(agriculture);
  
  // GEE images need a 'system:time_start' property to 
  // indicate the timestamp of the image
  // We copy it from the original image and return the result
  return classImage.copyProperties(image, ['system:time_start']);
};

// Apply the function on all images
var extractedCollection = landcoverFiltered.map(extractClasses);
print('Extracted Image Collection with Urban and Agriculture classes',
  extractedCollection);

// We modify the calculateNtl() function to use reduceRegions()
var calculateSol = function(image) {
  var year = image.date().get('year');
  var yearStart = ee.Date.fromYMD(year, 1, 1);
  var yearEnd = yearStart.advance(1, 'year');
  var classImage = extractedCollection
    .filter(ee.Filter.date(yearStart, yearEnd)).first();
  var urban = classImage.select('urban');
  var agriculture = classImage.select('agriculture');
  
  var urbanStatsCol = image.updateMask(urban).reduceRegions({
    collection: china,
    reducer: ee.Reducer.sum().setOutputs(['urban_sol']),
    scale: resolution,
    tileScale: 16
    })
  var agricultureStatsCol = image.updateMask(agriculture).reduceRegions({
    collection: urbanStatsCol,
    reducer: ee.Reducer.sum().setOutputs(['agriculture_sol']),
    scale: resolution,
    tileScale: 16})
  
  // ntl for the regions across all classes
  var statsCol = image.reduceRegions({
    collection: agricultureStatsCol,
    reducer: ee.Reducer.sum().setOutputs(['total_sol']),
    scale: resolution,
    tileScale: 16})

  // Process the results
  var results = statsCol.map(function(feature) {
    var solUrban = ee.List([feature.getNumber('urban_sol'), 0])
    .reduce(ee.Reducer.firstNonNull());
    var solAgriculture = ee.List([feature.getNumber('agriculture_sol'), 0])
    .reduce(ee.Reducer.firstNonNull());
    var solTotal = ee.List([feature.getNumber('total_sol'), 0])
    .reduce(ee.Reducer.firstNonNull());

    var province = feature.get('省');
    var province_id = feature.get('省代码');
    var city = feature.get('市');
    var city_id = feature.get('市代码');
    var county = feature.get('县');
    var county_id = feature.get('县代码');
    var newFeature = ee.Feature(null, {
      'urban_sol': solUrban,
      'agriculture_sol': solAgriculture,
      'total_sol': solTotal,
      'year': year,
      'province': province,
      'province_id': province_id,
      'city': city,
      'city_id': city_id,
      'county': county,
      'county_id': county_id
    });
    return newFeature;
  });
  
  return results;
};

// Apply this function on all images
// We get a nested collection since we have collection of collections
// Run flatten() to get a flat collection
var solByClassByRegionTimeSeries = filtered
  .map(calculateSol)
  .flatten();
print('Sum of Lights by Class by Region Time Series', solByClassByRegionTimeSeries);

// Export the results as a CSV

// Select the properties to export
var outputFields = ['province', 'province_id', 'city','city_id', 
                    'county','county_id', 'year',
                    'urban_sol', 'agriculture_sol', 'total_sol'];

Export.table.toDrive({
  collection: solByClassByRegionTimeSeries,
  description: 'SOL_Time_Series_by_county',
  folder: 'earthengine',
  fileNamePrefix: 'sol_by_year_by_class_by_county',
  fileFormat: 'CSV',
  selectors: outputFields});

For the day-time images

//Load china adm3 layer
var china = ee.FeatureCollection("projects/ee-chenyilin2025/assets/county");
var geometry = china.geometry();
Map.addLayer(china, {color: 'gray'}, 'Admin3 Boundaries');

//Load MODIS NPP and landcover images
var modisNpp = ee.ImageCollection('MODIS/061/MOD17A3HGF');
var band = ['Npp'];
var modisLandcover = ee.ImageCollection("MODIS/006/MCD12Q1");
// Using the IGBP classification
var classification = 'LC_Type1'

// NPP is measured in the unit 
// Kilogram Carbon per Square Meters per Year (kgC /m^2 / year)
// The MODIS images comes with a scale multiplier of 0.0001
// Multiplying the pixel values with scale gives us NPP
//
var scaleImage = function(image) {
  var scaled = image.multiply(0.0001)
  return scaled.copyProperties(image, ['system:index', 'system:time_start'])
};

// ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
//  Average Net Productivity by Class for all regions for all years
// ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
// Extract Landcover Classification for the year
var startYear = 2001;
var endYear = 2020;

var startDate = ee.Date.fromYMD(startYear, 1, 1);
var endDate = ee.Date.fromYMD(endYear + 1, 1, 1);

// Filter NPP image collection for the selected years 
// and select the band having class values
var nppFiltered = modisNpp
  .filter(ee.Filter.date(startDate, endDate))
  .filter(ee.Filter.bounds(geometry))
  .select(band);
print('Filtered Modis NPP Collection', nppFiltered);

// Apply the pixel scaling
var nppScaled = nppFiltered.map(scaleImage);
print('Filtered and Scaled Modis NPP Collection', nppScaled);

// Filter landcover image collection for the selected years 
// and select the band having class values
var landFiltered = modisLandcover
  .filter(ee.Filter.date(startDate, endDate))
  .select(classification);
print('Filtered MODIS Collection from ' + startYear + ' to ' + endYear,
  landFiltered);

var projection = nppFiltered.first().projection();
var resolution = projection.nominalScale();

// Write a function to extract urban and agriculture classes from
// a MODIS landcover image
var extractClasses = function(image) {
  var urban = image.eq(13)
    .rename('urban');
  var agriculture = image.eq(10)
    .or(image.eq(12))
    .or(image.eq(14))
    .rename('agriculture')
    
  // Create a 2-band image with the extracted classes
  var classImage = urban.addBands(agriculture);

  // GEE images need a 'system:time_start' property to 
  // indicate the timestamp of the image
  // We copy it from the original image and return the result
  return classImage.copyProperties(image, ['system:time_start']);
};

// Apply the function on all images
var extractedCollection = landFiltered.map(extractClasses);
print('Extracted Image Collection with Urban and Agriculture classes',
  extractedCollection);
  
// We now write a function to calculate ANPP for each class
// This function will return a ee.Feature() with ANPP as properties
var calculateANPP = function(image) {
  var year = image.date().get('year');
  var yearStart = ee.Date.fromYMD(year, 1, 1);
  var yearEnd = yearStart.advance(1, 'year');
  var classImage = extractedCollection
    .filter(ee.Filter.date(yearStart, yearEnd)).first();
  var urban = classImage.select('urban');
  var agriculture = classImage.select('agriculture');

var urbanStatsCol = image.updateMask(urban).reduceRegions({
    collection: china,
    reducer: ee.Reducer.mean().setOutputs(['urban_anpp']),
    scale: resolution,
    tileScale: 16
    })
    
  var agricultureStatsCol = image.updateMask(agriculture).reduceRegions({
    collection: urbanStatsCol,
    reducer: ee.Reducer.mean().setOutputs(['agriculture_anpp']),
    scale: resolution,
    tileScale: 16})

  // ANPP for the regions across all classes
  var statsCol = image.reduceRegions({
    collection: agricultureStatsCol,
    reducer: ee.Reducer.mean().setOutputs(['overall_anpp']),
    scale: resolution,
    tileScale: 16})

  // Process the results
  var results = statsCol.map(function(feature) {
    var anppUrban = ee.List([feature.getNumber('urban_anpp'), 0])
    .reduce(ee.Reducer.firstNonNull());
    var anppAgriculture = ee.List([feature.getNumber('agriculture_anpp'), 0])
    .reduce(ee.Reducer.firstNonNull());
    var anppOverall = ee.List([feature.getNumber('overall_anpp'), 0])
    .reduce(ee.Reducer.firstNonNull());

    var anppUrbanGrams = ee.Number(anppUrban).multiply(1000);
    var anppAgricultureGrams = ee.Number(anppAgriculture).multiply(1000);
    var anppOverallGrams = ee.Number(anppOverall).multiply(1000);

    var province = feature.get('省');
    var province_id = feature.get('省代码');
    var city = feature.get('市');
    var city_id = feature.get('市代码');
    var county = feature.get('县');
    var county_id = feature.get('县代码');
    var newFeature = ee.Feature(null, {
      'urban_anpp': anppUrbanGrams,
      'agriculture_anpp': anppAgricultureGrams,
      'overall_anpp': anppOverallGrams,
      'year': year,
      'province': province,
      'province_id': province_id,
      'city': city,
      'city_id': city_id,
      'county': county,
      'county_id': county_id
    });
    return newFeature;
  });
  
  return results;
};

// Apply this function on all images
// We get a nested collection since we have collection of collections
// Run flatten() to get a flat collection
var anppByClassByRegionTimeSeries = nppScaled
  .map(calculateANPP)
  .flatten();
print('ANPP by Class by Region Time Series', anppByClassByRegionTimeSeries);

// Export the results as a CSV

// Select the properties to export
var outputFields = ['province', 'province_id', 'city','city_id',
                    'year', 'county','county_id',
                    'urban_anpp', 'agriculture_anpp', 'overall_anpp'];

Export.table.toDrive({
  collection: anppByClassByRegionTimeSeries,
  description: 'ANPP_Time_Series_by_county',
  folder: 'earthengine',
  fileNamePrefix: 'anpp_by_year_by_class_by_county',
  fileFormat: 'CSV',
  selectors: outputFields});