User memory limit exceeded en Google Earth Engine: por qué ocurre y cómo resolverlo

/**
 * Experiment: what actually fixes "User memory limit exceeded"?
 *
 * Same calculation in five versions, from the typical naive script to the
 * corrected one, changing one thing at a time:
 *
 *     mean NDVI of 2024 (Sentinel-2 median, 10 m) for each municipality of
 *     Córdoba, Colombia, with reduceRegions.
 *
 * Each version runs twice:
 *     1. interactive (evaluate), where the error usually shows up: this
 *        script times each one and prints the result in the console;
 *     2. as an export task: click RUN on the five "memoria_..." tasks in the
 *        Tasks tab. When they finish, set READ_TASKS = true and run the script
 *        again to read their runtime and EECU-seconds.
 *
 * Keep this tab open until the console prints "Done". The interactive runs go
 * one after another and can take up to 25 minutes.
 */


/**
 * Settings
 */
var START = '2024-01-01';
var END = '2025-01-01';            // filterDate excludes END
var CRS = 'EPSG:32618';            // UTM 18N, the zone of Córdoba
var SCALE = 10;                    // Sentinel-2 B4 and B8
var DRIVE_FOLDER = 'earthengine';
var READ_TASKS = false;            // true: only read the finished exports


/**
 * Data
 *
 * GAUL 2025: GAUL 2015 is deprecated in the Earth Engine catalog and misses
 * the newer municipalities of Córdoba.
 */
var municipalities = ee.FeatureCollection('FAO/GAUL/2025/level2')
    .filter(ee.Filter.eq('ISO3_CODE', 'COL'))
    .filter(ee.Filter.inList('GAUL1_NAME', ['Córdoba', 'Cordoba']));

print('Polygons', municipalities.size());
print('Names', municipalities.aggregate_array('GAUL2_NAME').distinct().sort());

var region = municipalities.geometry();
Map.centerObject(municipalities, 8);
Map.addLayer(municipalities, {color: 'black'}, 'Municipalities');

var s2 = ee.ImageCollection('COPERNICUS/S2_SR_HARMONIZED')
    .filterBounds(region)
    .filterDate(START, END);

// SCL classes removed: 0 no data, 1 saturated, 3 cloud shadow,
// 8 and 9 cloud (medium and high probability), 10 cirrus.
var SCL_BAD = [0, 1, 3, 8, 9, 10];

function maskScl(image) {
  var bad = image.select('SCL').remap(SCL_BAD, SCL_BAD.map(function () { return 1; }), 0);
  return image.updateMask(bad.not());
}

// Median NDVI of the year. selectFirst keeps only B4, B8 and SCL before the
// median; otherwise the median runs over every band.
function ndviImage(selectFirst) {
  var collection = selectFirst ? s2.select(['B4', 'B8', 'SCL']) : s2;
  return collection.map(maskScl).median()
      .normalizedDifference(['B8', 'B4'])
      .rename('NDVI');
}


/**
 * Versions: each one changes a single thing from the previous one
 */
var VERSIONS = [
  {
    key: 'v1_ingenuo',
    label: 'Todas las bandas, clip() y reproject(), tileScale 1',
    // What a typical script (or an AI assistant) writes: every band,
    // clip() to the region and reproject() "to fix the scale".
    image: function () { return ndviImage(false).clip(region).reproject({crs: CRS, scale: SCALE}); },
    tileScale: 1
  },
  {
    key: 'v2_sin_clip',
    label: 'Sin clip() ni reproject()',
    image: function () { return ndviImage(false); },
    tileScale: 1
  },
  {
    key: 'v3_bandas',
    label: 'Solo B4, B8 y SCL antes de la mediana',
    image: function () { return ndviImage(true); },
    tileScale: 1
  },
  {
    key: 'v4_tilescale4',
    label: 'Lo anterior y tileScale 4',
    image: function () { return ndviImage(true); },
    tileScale: 4
  },
  {
    key: 'v5_tilescale16',
    label: 'Lo anterior y tileScale 16',
    image: function () { return ndviImage(true); },
    tileScale: 16
  }
];

function perMunicipality(version) {
  return version.image().reduceRegions({
    collection: municipalities,
    reducer: ee.Reducer.mean(),
    crs: CRS,
    scale: SCALE,
    tileScale: version.tileScale
  }).select(['GAUL2_NAME', 'mean']);
}


/**
 * 1. Exports: they appear in the Tasks tab. Click RUN on each one.
 */
VERSIONS.forEach(function (version) {
  if (READ_TASKS) return;
  Export.table.toDrive({
    collection: perMunicipality(version),
    description: 'memoria_' + version.key,
    folder: DRIVE_FOLDER,
    fileFormat: 'CSV'
  });
});


/**
 * 2. Interactive runs, one at a time
 *
 * From the corrected version to the naive one: if Earth Engine reuses any
 * cached piece, the advantage goes to the naive versions, not to the fixes.
 */
var order = VERSIONS.slice().reverse();
var rows = [];
var values = {};    // key -> {municipality: mean NDVI}

function runNext(i) {
  if (i >= order.length) {
    showTables();
    return;
  }
  var version = order[i];
  var t0 = Date.now();
  print('Running ' + version.key + '…');
  perMunicipality(version).evaluate(function (result, error) {
    var seconds = Math.round((Date.now() - t0) / 100) / 10;
    var state = error ? 'FAILED' : 'COMPLETED';
    if (!error) {
      values[version.key] = {};
      result.features.forEach(function (f) {
        values[version.key][f.properties.GAUL2_NAME] = f.properties.mean;
      });
    }
    print(version.key + ': ' + state + ', ' + seconds + ' s' + (error ? ' · ' + error : ''));
    rows.push([version.key, version.label, state, seconds, error || '']);
    runNext(i + 1);
  });
}


/**
 * 3. Results as tables: open each one with the arrow (↗) at its top right
 *    and download it as CSV.
 */
function showTables() {
  // Interactive runs, in the order of VERSIONS
  rows.sort(function (a, b) { return a[0] < b[0] ? -1 : 1; });
  var runs = {
    cols: [
      {label: 'version', type: 'string'},
      {label: 'cambio', type: 'string'},
      {label: 'estado', type: 'string'},
      {label: 'segundos', type: 'number'},
      {label: 'error', type: 'string'}
    ],
    rows: rows.map(function (r) {
      return {c: r.map(function (v) { return {v: v}; })};
    })
  };
  print(ui.Chart(runs, 'Table').setOptions({title: 'experimento_memoria_interactivo'}));

  // Mean NDVI per municipality of every interactive run that finished
  var keys = VERSIONS.map(function (v) { return v.key; })
      .filter(function (k) { return values[k]; });
  if (keys.length) {
    var names = Object.keys(values[keys[0]]).sort();
    var ndvi = {
      cols: [{label: 'municipio', type: 'string'}].concat(keys.map(function (k) {
        return {label: k, type: 'number'};
      })),
      rows: names.map(function (m) {
        return {c: [{v: m}].concat(keys.map(function (k) {
          return {v: values[k][m]};
        }))};
      })
    };
    print(ui.Chart(ndvi, 'Table').setOptions({title: 'ndvi_por_municipio'}));
  }
  print('Done');
}


/**
 * 4. Runtime and EECU-seconds of the exports (READ_TASKS = true)
 *
 * Run it when the five tasks have finished. If an export was launched more
 * than once, it takes the latest run.
 */
function readTasks() {
  ee.data.listOperations(100, function (ops, error) {
    if (error) {
      print('Error: ' + error);
      return;
    }

    // Latest run of each task, in case an export was launched more than once
    var latest = {};
    ops.forEach(function (op) {
      var m = op.metadata || {};
      if (!m.description || m.description.indexOf('memoria_') !== 0) return;
      if (!latest[m.description] || m.createTime > latest[m.description].metadata.createTime) {
        latest[m.description] = op;
      }
    });

    var rows = Object.keys(latest).sort().map(function (name) {
      var op = latest[name];
      var m = op.metadata;
      var end = m.endTime || m.updateTime;
      var seconds = m.startTime && end ?
          Math.round((Date.parse(end) - Date.parse(m.startTime)) / 100) / 10 : null;
      var error = op.error ? op.error.message : '';
      return {c: [
        {v: name.replace('memoria_', '')},
        {v: m.state},
        {v: seconds},
        {v: m.batchEecuUsageSeconds !== undefined ? Number(m.batchEecuUsageSeconds) : null},
        {v: error}
      ]};
    });

    if (!rows.length) {
      print('No "memoria_..." tasks found.');
      return;
    }

    print(ui.Chart({
      cols: [
        {label: 'version', type: 'string'},
        {label: 'estado', type: 'string'},
        {label: 'segundos', type: 'number'},
        {label: 'eecu_segundos', type: 'number'},
        {label: 'error', type: 'string'}
      ],
      rows: rows
    }, 'Table').setOptions({title: 'experimento_memoria_exportacion'}));
  });
}


if (READ_TASKS) {
  readTasks();
} else {
  runNext(0);
}

"""
Experiment: what actually fixes "User memory limit exceeded"?

Same calculation in five versions, from the typical naive script to the
corrected one, changing one thing at a time:

    mean NDVI of 2024 (Sentinel-2 median, 10 m) for each municipality of
    Córdoba, Colombia, with reduceRegions.

Each version runs twice:
    1. interactive (getInfo), where the error usually shows up;
    2. as an export task, which reports its runtime and EECU-seconds.

Output: experimento_memoria.csv with one row per version and mode
(state, error message, seconds, EECU-seconds), plus the per-municipality
NDVI of every run that finished, to check that the fixes do not change the
numbers.

Run it in Colab or Jupyter. It takes between 20 and 60 minutes.
"""

import csv
import time

import ee

ee.Authenticate()
ee.Initialize(project='tu-proyecto-gcp')


# ---------------------------------------------------------------------------
# Settings
# ---------------------------------------------------------------------------

START, END = '2024-01-01', '2025-01-01'   # filterDate excludes END
CRS = 'EPSG:32618'                         # UTM 18N, the zone of Córdoba
SCALE = 10                                 # Sentinel-2 B4 and B8
DRIVE_FOLDER = 'earthengine'
OUTPUT_CSV = 'experimento_memoria.csv'
POLL_SECONDS = 30


# ---------------------------------------------------------------------------
# Data
# ---------------------------------------------------------------------------

# GAUL 2025: GAUL 2015 is deprecated in the Earth Engine catalog and misses
# the newer municipalities of Córdoba.
municipalities = (
    ee.FeatureCollection('FAO/GAUL/2025/level2')
    .filter(ee.Filter.eq('ISO3_CODE', 'COL'))
    .filter(ee.Filter.inList('GAUL1_NAME', ['Córdoba', 'Cordoba']))
)
print('Polygons:', municipalities.size().getInfo(),
      sorted(set(municipalities.aggregate_array('GAUL2_NAME').getInfo())))
region = municipalities.geometry()

s2 = (
    ee.ImageCollection('COPERNICUS/S2_SR_HARMONIZED')
    .filterBounds(region)
    .filterDate(START, END)
)

# SCL classes removed: 0 no data, 1 saturated, 3 cloud shadow,
# 8 and 9 cloud (medium and high probability), 10 cirrus.
SCL_BAD = [0, 1, 3, 8, 9, 10]


def mask_scl(image):
    scl = image.select('SCL')
    bad = scl.remap(SCL_BAD, [1] * len(SCL_BAD), 0)
    return image.updateMask(bad.Not())


def ndvi_image(select_first):
    """Median NDVI of the year. select_first keeps only B4, B8 and SCL
    before the median; otherwise the median runs over every band."""
    collection = s2.select(['B4', 'B8', 'SCL']) if select_first else s2
    median = collection.map(mask_scl).median()
    return median.normalizedDifference(['B8', 'B4']).rename('NDVI')


# ---------------------------------------------------------------------------
# Versions: each one changes a single thing from the previous one
# ---------------------------------------------------------------------------

def v1_naive():
    # What a typical script (or an AI assistant) writes: every band,
    # clip() to the region and reproject() "to fix the scale".
    return ndvi_image(False).clip(region).reproject(crs=CRS, scale=SCALE), 1


def v2_no_clip_reproject():
    return ndvi_image(False), 1


def v3_select_bands():
    return ndvi_image(True), 1


def v4_tilescale_4():
    return ndvi_image(True), 4


def v5_tilescale_16():
    return ndvi_image(True), 16


VERSIONS = [
    ('v1_ingenuo', 'Todas las bandas, clip() y reproject(), tileScale 1', v1_naive),
    ('v2_sin_clip', 'Sin clip() ni reproject()', v2_no_clip_reproject),
    ('v3_bandas', 'Solo B4, B8 y SCL antes de la mediana', v3_select_bands),
    ('v4_tilescale4', 'Lo anterior y tileScale 4', v4_tilescale_4),
    ('v5_tilescale16', 'Lo anterior y tileScale 16', v5_tilescale_16),
]


def per_municipality(version_fn):
    image, tile_scale = version_fn()
    return image.reduceRegions(
        collection=municipalities,
        reducer=ee.Reducer.mean(),
        crs=CRS,
        scale=SCALE,
        tileScale=tile_scale,
    ).select(['GAUL2_NAME', 'mean'])


# ---------------------------------------------------------------------------
# 1. Exports: start all of them now, they run in the background
# ---------------------------------------------------------------------------

tasks = {}
for key, _, fn in VERSIONS:
    task = ee.batch.Export.table.toDrive(
        collection=per_municipality(fn),
        description=f'memoria_{key}',
        folder=DRIVE_FOLDER,
        fileFormat='CSV',
    )
    task.start()
    tasks[key] = task
    print('Export started:', key)


# ---------------------------------------------------------------------------
# 2. Interactive runs, one at a time
#
# From the corrected version to the naive one: if Earth Engine reuses any
# cached piece, the advantage goes to the naive versions, not to the fixes.
# ---------------------------------------------------------------------------

rows = []
values = {}   # key -> {municipality: mean NDVI}

for key, label, fn in reversed(VERSIONS):
    print('Interactive:', key, '...', end=' ', flush=True)
    t0 = time.time()
    try:
        result = per_municipality(fn).getInfo()
        state, error = 'COMPLETED', ''
        values[key] = {
            f['properties']['GAUL2_NAME']: f['properties'].get('mean')
            for f in result['features']
        }
    except ee.EEException as e:
        state, error = 'FAILED', str(e)
    seconds = round(time.time() - t0, 1)
    print(state, seconds, 's', error[:80])
    rows.append({
        'version': key, 'cambio': label, 'modo': 'interactivo',
        'estado': state, 'error': error, 'segundos': seconds,
        'eecu_segundos': '',
    })


# ---------------------------------------------------------------------------
# 3. Wait for the exports and read their runtime and EECU usage
# ---------------------------------------------------------------------------

pending = set(tasks)
while pending:
    for key in sorted(pending):
        if tasks[key].status()['state'] in ('COMPLETED', 'FAILED', 'CANCELLED'):
            pending.discard(key)
    if pending:
        print('Exports still running:', ', '.join(sorted(pending)))
        time.sleep(POLL_SECONDS)

for key, label, _ in VERSIONS:
    status = tasks[key].status()
    start = status.get('start_timestamp_ms')
    end = status.get('update_timestamp_ms')
    seconds = round((end - start) / 1000, 1) if start and end else ''
    rows.append({
        'version': key, 'cambio': label, 'modo': 'exportacion',
        'estado': status['state'], 'error': status.get('error_message', ''),
        'segundos': seconds,
        'eecu_segundos': status.get('batch_eecu_usage_seconds', ''),
    })


# ---------------------------------------------------------------------------
# 4. Do the fixes change the numbers? Largest difference in mean NDVI
#    against the first version that finished interactively.
# ---------------------------------------------------------------------------

reference_key = next((k for k, _, _ in VERSIONS if k in values), None)
for row in rows:
    key = row['version']
    if row['modo'] == 'interactivo' and key in values and reference_key:
        ref = values[reference_key]
        diffs = [
            abs(v - ref[m]) for m, v in values[key].items()
            if v is not None and ref.get(m) is not None
        ]
        row['max_dif_ndvi_vs_' + reference_key] = round(max(diffs), 5) if diffs else ''


# ---------------------------------------------------------------------------
# Save
# ---------------------------------------------------------------------------

fields = ['version', 'cambio', 'modo', 'estado', 'segundos', 'eecu_segundos', 'error']
fields += sorted({f for r in rows for f in r} - set(fields))
with open(OUTPUT_CSV, 'w', newline='', encoding='utf-8') as f:
    writer = csv.DictWriter(f, fieldnames=fields)
    writer.writeheader()
    writer.writerows(rows)

with open('ndvi_por_municipio.csv', 'w', newline='', encoding='utf-8') as f:
    keys = [k for k, _, _ in VERSIONS if k in values]
    writer = csv.writer(f)
    writer.writerow(['municipio'] + keys)
    for m in sorted({m for k in keys for m in values[k]}):
        writer.writerow([m] + [values[k].get(m) for k in keys])

print('\nSaved', OUTPUT_CSV, 'and ndvi_por_municipio.csv')

# In Colab the local disk is deleted when the session ends: download both files.
try:
    from google.colab import files
    files.download(OUTPUT_CSV)
    files.download('ndvi_por_municipio.csv')
except ImportError:
    pass
for r in rows:
    print(f"{r['version']:16} {r['modo']:12} {r['estado']:10} "
          f"{r['segundos']!s:>8} s  {r['eecu_segundos']!s:>8} EECU-s  {r['error'][:60]}")