Is the lake or river my town drinks from getting choked with toxic algae?

Cyanobacteria (you may know them as “blue-green algae”) are tiny organisms that can explode in number in warm, calm fresh water, turning a lake green and scummy. When they bloom in the lake or reservoir a town drinks from, they can poison the water supply and shut down fishing. A satellite can spot the green tint of a bloom from orbit, so you can watch one build up over your specific waterbody even when nobody is out there testing it.

Which data to use, and how

Use this dataset: CyAN cyanobacteria index (short name MERGED_S3_OLCI_L3m_CYAN). It takes pictures from the Sentinel-3 satellite and turns them into a ready-made “how much blue-green algae is here” score for every spot on a lake, covering 2016 → today. The catch: it’s pre-computed for US lakes only. If your lake is outside the US, you’d instead start from the raw scenes (OLCI inland-water Level-2, OLCIS3A_L2_ILW / OLCIS3B_L2_ILW) and compute the score yourself. That’s harder, so begin with the ready-made CyAN product.

Then do four steps (the runnable code further down does all of this; this is just the idea):

1. Pick your lake. Draw a small latitude/longitude box around it, and if you want, mark the exact spot where the drinking-water intake pipe sits.
2. Read the algae score for each picture. The data comes as one image roughly every 7 days. Each pixel holds a number that means “how much blue-green algae is here” (more on that number below).
3. Boil each image down to a few numbers: the average score over the lake, how big the bloom is (in km²), and the score right at the intake pipe.
4. Watch it over time. Build a “normal” seasonal pattern from past years so you know your usual bloom window, then flag any date when a strong bloom sits near the intake. That’s the moment a water utility most needs the warning.

That’s the whole method. Below, “How a scientist answers this” names the exact tools, and the Run it cell does it live on synthetic data so you can see each step work.

What you can find out

• Is a bloom forming on your source water right now? The index maps blue-green algae pixel-by-pixel across the whole lake or reservoir.
• How bad and how big it is. Sort the score into low / moderate / high, and measure the bloom’s area in km².
• Is it near the drinking-water intake? Read the score right at the intake pipe and the pixels around it, not just the lake-wide average.
• When it normally peaks. Stack up the years (2016+) to learn your usual bloom window, so an early or unusually fierce season stands out.
• Whether it’s getting worse year to year. Track the peak bloom size and intensity across the whole Sentinel-3 record for a trend.

What it can’t tell you

• Lakes outside the US (with the easy dataset). The ready-made CyAN index only covers US lakes and reservoirs. Elsewhere you have to start from the raw OLCI inland-water scenes and build the index yourself.
• The toxin level. The index measures how many blue-green algae are present, not how much poison (microcystin, cylindrospermopsin) they’re making. To get toxin amounts, someone has to scoop a water sample and run a lab test.
• Whether your tap water is safe. The satellite sees the raw lake water, not what comes out after the treatment plant has cleaned it. You’d pair this with the utility’s own tap-water testing.
• Small ponds, narrow rivers, or the near-shore strip. Each pixel is ~300 m across, and pixels touching land are thrown out, so many small reservoirs are too small to see.
• Blooms hiding below the surface. A satellite only sees the thin top layer of water. A bloom mixed deep down or sitting on the bottom is invisible to it.
• Which species it is. The index flags a “looks like blue-green algae” signal, not the exact species or whether it’s a toxin-making strain. Only a lab can confirm that.

Gotchas to watch for

• The pixel value is a code, not a toxin amount. Each pixel is a number from 0 to 255: 0 means “nothing detected”, 1–253 is the algae index climbing higher, 254 means “this is land”, and 255 means “no data here”. So before you do math, keep only the 1–253 water pixels and treat the rest as blanks. Otherwise “land” (254) would inflate your averages.
• Index is not the same as toxin. A high index means “go test the water now,” not “the poison level is X.” Always back up a satellite alert with an actual water sample.
• Clouds and wind. A cloudy sky hides the lake completely, and wind can shove a whole bloom against one shoreline within hours. Treat each clear image as a single snapshot, not a continuous video.
• It’s coarse. At ~300 m per pixel, with shoreline pixels removed, small reservoirs and skinny river stretches come out blurry or missing entirely.
• Muddy or stained water can fool it (Case-II water). Inland water full of sediment or dissolved plant matter (“Case-II” water, the optically messy kind) can confuse the reading. The index is tuned for blue-green algae but isn’t perfect, which is why a lab check matters.

The real-data code

The Run it cell above runs this method on synthetic data with no login. Below is the whole recipe against the real archive: search the CyAN composites, open each GeoTIFF, score the bloom over your lake and at the intake, build a seasonal climatology, print a verdict, and plot it. It needs a free Earthdata Login. The merged CyAN product ships as Cloud-Optimized GeoTIFFs (image files, one 7-day composite each), gridded to a US-wide map projection (EPSG:5070) at 300 m. So it’s US lakes only, and you read it with rioxarray, not plain xarray.

import earthaccess
import rioxarray
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from rasterio.warp import transform_bounds

earthaccess.login(strategy="netrc")

# Western Lake Erie — Toledo drinking-water intake region (lon/lat box)
aoi = (-83.5, 41.4, -82.5, 42.0)
intake_lon, intake_lat = -83.26, 41.69   # approx. Toledo Collins Park WTP crib intake

# 1. Merged Sentinel-3 OLCI cyanobacteria index (CyAN), 7-day composites, 2016+
results = earthaccess.search_data(
    short_name="MERGED_S3_OLCI_L3m_CYAN",
    bounding_box=aoi,
    temporal=("2016-04-01", "2025-10-31"),   # full record; bloom season is summer
)
files = earthaccess.open(results)

# 2. Open each GeoTIFF and boil it down to a bloom-severity record for the lake AOI.
records, last_map = [], None
for f in files:
    da = rioxarray.open_rasterio(f, masked=False).squeeze()      # DN raster, EPSG:5070
    # reproject the lon/lat box into the image's projection, then crop to it
    xmin, ymin, xmax, ymax = transform_bounds("EPSG:4326", da.rio.crs, *aoi)
    box = da.rio.clip_box(xmin, ymin, xmax, ymax)
    win = box.values.astype("float32")

    cyano = (win >= 1) & (win <= 253)             # real water pixels with a reading
    dn = np.where(cyano, win, np.nan)             # algae value where detectable, else blank

    # value at the intake pixel (reproject just the crop back to lon/lat to index it)
    iv = box.rio.reproject("EPSG:4326").sel(
        x=intake_lon, y=intake_lat, method="nearest").values
    intake_dn = float(iv) if 1 <= iv <= 253 else np.nan

    records.append({
        "date":          pd.to_datetime(f.granule.get("time_start", None), errors="coerce"),
        "mean_dn":       float(np.nanmean(dn)) if cyano.any() else 0.0,
        "max_dn":        float(np.nanmax(dn))  if cyano.any() else 0.0,
        "bloom_area_px": int(cyano.sum()),       # any detectable algae
        "high_area_px":  int((win >= 200).sum()),# high end of the index
        "intake_dn":     intake_dn,
    })
    last_map = dn                                  # keep the latest scene for the map panel

ts = pd.DataFrame(records).dropna(subset=["date"]).sort_values("date")
ts["bloom_area_km2"] = ts["bloom_area_px"] * (0.3 * 0.3)   # ~0.09 km² per 300 m pixel

# 3. Seasonal climatology + intake alerting (a high index at the intake → flag the date).
ts["intake_alert"] = ts["intake_dn"] >= 200
climatology = ts.groupby(ts["date"].dt.month)["bloom_area_km2"].mean()
peak_month = int(climatology.idxmax())
worst = ts.loc[ts["bloom_area_km2"].idxmax()]
n_alerts = int(ts["intake_alert"].sum())

# 4. Verdict + plots: lake-wide bloom area over time, the seasonal cycle, and the latest map.
print(f"Worst bloom: {worst['bloom_area_km2']:.0f} km^2 on {worst['date'].date()}; "
      f"typical peak in month {peak_month}; {n_alerts} dates with a strong bloom at the intake.")

fig, ax = plt.subplots(1, 3, figsize=(15, 4))
ax[0].plot(ts["date"], ts["bloom_area_km2"]); ax[0].set_title("Bloom area (km²)")
ax[0].scatter(ts.loc[ts["intake_alert"], "date"],
              ts.loc[ts["intake_alert"], "bloom_area_km2"], c="r", label="intake alert")
ax[0].legend()
climatology.plot.bar(ax=ax[1], title="Mean bloom area by month")
im = ax[2].imshow(last_map, vmin=1, vmax=253, cmap="viridis")
ax[2].set_title("Latest CyAN index"); fig.colorbar(im, ax=ax[2])
plt.tight_layout(); plt.show()

# Before you trust it: a high index means "go test the water now," not a toxin level — confirm any
# alert with an actual water sample (microcystin lab test), per the gotchas above.

What you’ll get out of this: a bloom-severity time series for your source water (lake-wide bloom area in km², plus the index right at the intake), a seasonal climatology of the normal bloom window (for western Lake Erie, late July–September), an intake alert flag for the dates a strong bloom sat near the intake, and a map of the bloom for any flagged date.

Where the data comes from

All from NASA through one free login (Earthdata Login). The ready-made CyAN index and the raw OLCI inland-water scenes both come from NASA’s Ocean Biology distribution center (OB.DAAC), served through the CyAN project (Cyanobacteria Assessment Network), which exists specifically to turn satellite color into algae warnings for inland US drinking-water lakes.

Sources

• CyAN (Cyanobacteria Assessment Network) project: https://oceancolor.gsfc.nasa.gov/about/projects/cyan/
• CyAN data on OB.DAAC OceanColor: https://oceancolor.gsfc.nasa.gov/
• Sentinel-3 OLCI mission: https://sentinels.copernicus.eu/web/sentinel/missions/sentinel-3
• EPA CyAN drinking-water / public-health use: https://www.epa.gov/water-research/cyanobacteria-assessment-network-cyan