geeViz.fsInsights - FIA and LCMS, made usable

Two public USDA Forest Service APIs that answer complementary halves of the same question and look nothing alike:

LCMS

FIA

Nature

Wall-to-wall 30 m map

Probability sample of plots

Span

1985-2025, annual

1984-present, multi-year panels

Answers

What changed, and where

What is there, and how much

Uncertainty

Map accuracy

Design-based standard error

Neither needs authentication. Earth Engine is optional - only the arbitrary-geometry path uses it.

Two things to know before trusting a number from here:

  1. Every FIA estimate carries its sampling error, and cells too thin to report are flagged. An FIA estimate without its error is not a fact.

  2. LCMS and FIA are not directly comparable. This package makes the comparison easy and labels it; it will not hand you a blended number.

github github

1. Discovery works offline

FIA’s /fullreport takes an attribute, two groupings, and an evaluation: 752 x 96 x 96 x 1129. Nobody memorizes that, and the docs still say ‘under construction’.

All three catalogs ship in the wheel, so search is instant and needs no network. That is not theoretical - FIADB-API returned HTTP 200 with an HTML error page across every endpoint while this package was being written, and these cells kept working.

from geeViz import fsInsights as fs

# 752 estimate attributes - search them
fs.find_attributes('carbon', limit=8)
# Narrow by land basis and the evaluation type an attribute requires
fs.find_attributes('net growth volume', land_basis='Timberland', limit=6)
# 1,129 evaluations. growth_only=True keeps those supporting growth,
# removals and mortality - a prerequisite the docs do not surface.
fs.find_evaluations('Oregon')
# 96 grouping variables, with the real database column behind each
fs.find_groupings('species', limit=6)

describe_grouping returns FIA’s own prose for a variable, including what each code value means - the difference between a column of integers and a column you can interpret.

print(fs.describe_grouping('Aspect')[:600])

2. FIA estimates, with their error attached

This is what the package is built around. Every cell comes back with se_pct and plots, and is flagged when SE% > 30 or n < 30.

n = 0 and n = 1 get their own reason codes rather than folding into ‘high error’, because a zero-plot cell reports a zero standard error - which reads as maximum precision when it means no information.

# Local validation first: this pairing is impossible and never leaves.
try:
    fs.validate(wc=12025, snum=999999)
except fs.FIAValidationError as e:
    print('caught locally:', e)

A note on FIADB-API availability

The estimate endpoint is intermittently unavailable, and it fails in an unusual way: HTTP 200 with an HTML error page rather than a 5xx (‘Internal Server Error: list index out of range’). During development the same query succeeded and failed minutes apart.

fsInsights detects that, retries it like the 5xx it is, and raises UpstreamError with the real message. The next cell probes once so the rest of the notebook can degrade gracefully instead of erroring out - everything above this point is offline and unaffected.

FIA_UP = True
try:
    _probe = fs.estimate(wc=12025, snum=2)
    print('FIADB-API is up -', len(_probe), 'cells from the probe')
except fs.UpstreamError as e:
    FIA_UP = False
    print('FIADB-API is down right now:')
    print(' ', str(e)[:180])
    print('\nEstimate cells below will skip. Discovery still works.')
# Forest area by county x forest type group, Alabama
df = None
if FIA_UP:
    df = fs.estimate(wc=12025, snum=2,
                     rselected='County code and name',
                     cselected='Forest type group')
    print(f'{len(df)} cells, units = ' + str(df['units'].iloc[0]))
    cols = ['column', 'estimate', 'se_pct', 'plots', 'unreliable']
    display(df[df['row'] == 'Total'][cols].head(8))
else:
    print('skipped - FIADB-API unavailable')

Look at the flagged rows. On a real cross-tabulation like this one, a large fraction of cells are too thin to report - the honest shape of a sample split many ways, and worth seeing rather than smoothing over.

When this was last run against a live API: 380 of 499 cells flagged.

if df is not None:
    flagged = df[df['unreliable']]
    print(f'{len(flagged)} of {len(df)} cells flagged unreliable')
    display(flagged[['column', 'estimate', 'se_pct', 'plots',
                     'unreliable_reason']].head(5))
else:
    print('skipped - no estimate to inspect')
# reliable() drops them - a caller decision, never done silently
if df is not None:
    clean = fs.reliable(df)
    print(f'{len(clean)} cells survive the default thresholds')
else:
    print('skipped')

3. LCMS - land cover, land use, and change

Small, clean, release-versioned. 3,643 precomputed summary areas: 3,137 counties, 502 ranger districts, and CONUS / All-Lands rollups.

Releases are not interchangeable

Two traps worth seeing before you pin a release:

  1. Products differ. 2025-6 is a tree canopy release carrying only NLCD_Percent_Tree_Canopy_Cover - no Land_Cover at all. So ‘the latest release’ is ambiguous unless you say latest of what.

  2. Study-area coverage is not monotonic. 2024-10 covers CONUS, AK, HAWAII and PRUSVI; the newer 2025-11 covers only CONUS and AK. Work in Hawaii or Puerto Rico has to pin an older release.

for r in fs.lcms_releases():
    prods = [p['Name'] for p in r['Products']]
    sas = [s['StudyArea'] for s in r['StudyAreas']]
    print(f"{r['VersionNumber']:9s} {r['StartYear']}-{r['EndYear']} "
          f"areas={r['SummaryAreaCount']:5d}")
    print(f"          products={prods}")
    print(f"          study_areas={sas}")
# 'latest' resolves per product
print('latest with Land_Cover :', fs.latest_release('Land_Cover'))
print('latest with TCC        :',
      fs.latest_release('NLCD_Percent_Tree_Canopy_Cover'))
# Asking a tree-canopy release for land cover is caught locally,
# and the error names a release that would work.
try:
    fs.lcms_summary('Land_Cover', state='Oregon', county='Crook',
                    year=2024, release='2025-6')
except ValueError as e:
    print(e)
# Class areas for one county, one year
lc = fs.lcms_summary(state='Oregon', county='Crook', year=2024)
lc[['year', 'class_name', 'acres', 'source']]
# Omit `year` for the full 1985-2025 series (~64 KB - fetch eagerly)
ts = fs.lcms_summary('Change', state='Oregon', county='Crook')
print(f'{len(ts)} rows across ' + str(ts['year'].nunique()) + ' years')
ts[ts['class_name'] == 'Wildfire'].nlargest(5, 'acres')[
    ['year', 'class_name', 'acres']]

Palettes come from the API

Each class arrives with its ClassValue and ClassPalette, so a legend cannot drift out of sync with a release. Any hardcoded LCMS palette is a maintenance liability this removes.

fs.lcms_classes('Change').head(8)
vp = fs.lcms_vis_params('Change')
print('min=', vp['min'], 'max=', vp['max'],
      'palette entries=', len(vp['palette']))
list(vp['classLegendDict'].items())[:5]

4. Arbitrary geometry (needs Earth Engine)

The API serves only its 3,643 precomputed areas - bbox, geojson and geometry are all rejected and POST returns 403. So lcms_summary dispatches a geometry to Earth Engine instead, computing the same class areas over any polygon.

Same call, same output shape, and the frame records which backend produced it. Skip this section if Earth Engine is not initialized - everything above works without it.

RUN_EE = False  # flip to True once ee is initialized

if RUN_EE:
    import ee
    import geeViz.geeView as gv
    aoi = ee.Geometry.Rectangle([-120.9, 44.1, -120.5, 44.4])
    ee_lc = fs.lcms_summary('Land_Cover', geometry=aoi,
                            year=2024, scale=90)
    display(ee_lc[['year', 'class_name', 'acres', 'source']])
else:
    print('RUN_EE is False - skipping the Earth Engine path')

For display of LCMS in the geeViz viewer, prefer the built-in thematic handling: the assets carry class names, values and palettes, so autoViz reads them directly.

if RUN_EE:
    from geeViz.fsInsights.lcms_ee import lcms_ee_collection
    Map = gv.Map
    Map.clearMap()
    lc2024 = lcms_ee_collection('Land_Cover').filter(
        ee.Filter.eq('year', 2024))
    Map.addLayer(lc2024.mosaic(), {'autoViz': True},
                 'LCMS Land Cover 2024')
    Map.centerObject(aoi)
    Map.view()
else:
    print('RUN_EE is False - skipping the map')

5. Putting them side by side - carefully

align.compare_area reports both figures and their difference, and returns the caveats with the numbers.

The difference is an observation, not an error term. A 15% gap does not mean either source is 15% wrong: FIA’s ‘forest land’ is a land-use definition based on stocking and potential, while LCMS land cover describes present canopy - so a recently harvested stand stays forest land in FIA while LCMS may map it as grass the same year.

from geeViz.fsInsights import align

result = None
if FIA_UP:
    result = align.compare_area(wc=412022, state='Oregon',
                                county='Crook', year=2024)
    print(align.summarize_comparison(result))
else:
    print('skipped - the FIA half of the comparison is unavailable')

# The LCMS half stands on its own either way:
treed = align.lcms_tree_area(state='Oregon', county='Crook', year=2024)
treed[['year', 'acres', 'source', 'estimator']]

The caveats travel with the numbers. A number this easy to quote is a number that travels without its footnotes unless you attach them.

# Available with or without a live FIA call - these are static
# facts about the two datasets, not properties of one comparison.
for c in align.COMPARISON_CAVEATS:
    print('-', c)
    print()

6. Cache management

Vocabularies are cached under ~/.geeViz/fsInsights with a 30-day stale-while-revalidate, falling back to the bundled snapshot and then to an empty result. A lookup never blocks on the network and never raises.

refresh_all() is the escape hatch for ‘a new release just dropped’ - not the normal path.

print('cache dir:', fs.cache_dir())
# fs.refresh_all()   # uncomment to force a refresh