Download notebook Run in Google Colab

Multi-band feature extraction¶

light-curve supports multi-band light curves natively in every feature extractor. Pass bands=["g", "r", ...] to the constructor to switch a feature into multiband mode: the feature is evaluated independently per passband and the outputs are concatenated, with names like amplitude_g, amplitude_r.

This tutorial covers:

1. Building a multiband light curve
2. Per-band feature extraction
3. Pure-multiband color features (ColorOfMedian, ColorOfMaximum, ColorOfMinimum, ColorSpread)
4. Mixed single-band + multiband Extractor
5. Multiband Bins and Periodogram
6. Batch processing with .many()
7. Arrow input for multiband features

In [1]:

# %pip install light-curve pyarrow

# %pip install light-curve pyarrow

In [2]:

import light_curve as licu
import numpy as np

rng = np.random.default_rng(42)

import light_curve as licu import numpy as np rng = np.random.default_rng(42)

1. Building a multiband light curve¶

A multiband light curve is just a regular (t, m, sigma) triple plus a band string array that labels each observation. We sort everything by time so we can pass sorted=True later.

In [3]:

def make_multiband_lc(band_labels, n_per_band=80, rng=None):
    """Return time-sorted (t, m, sigma, band) with n_per_band obs per passband."""
    rng = np.random.default_rng(rng)
    parts = []
    for label in band_labels:
        t_b = rng.uniform(0, 100, n_per_band)
        m_b = rng.normal(15.0, 0.3, n_per_band)
        s_b = np.full(n_per_band, 0.05)
        b_b = np.full(n_per_band, label)
        parts.append((t_b, m_b, s_b, b_b))
    t, m, s, b = (np.concatenate([p[i] for p in parts]) for i in range(4))
    idx = np.argsort(t)
    return t[idx], m[idx], s[idx], b[idx]


band_labels = ["g", "r", "i"]
t, m, sigma, band = make_multiband_lc(band_labels, rng=rng)
print(f"Total observations: {len(t)}")
print(f"Band counts: { {b: int(np.sum(band == b)) for b in band_labels} }")

def make_multiband_lc(band_labels, n_per_band=80, rng=None): """Return time-sorted (t, m, sigma, band) with n_per_band obs per passband.""" rng = np.random.default_rng(rng) parts = [] for label in band_labels: t_b = rng.uniform(0, 100, n_per_band) m_b = rng.normal(15.0, 0.3, n_per_band) s_b = np.full(n_per_band, 0.05) b_b = np.full(n_per_band, label) parts.append((t_b, m_b, s_b, b_b)) t, m, s, b = (np.concatenate([p[i] for p in parts]) for i in range(4)) idx = np.argsort(t) return t[idx], m[idx], s[idx], b[idx] band_labels = ["g", "r", "i"] t, m, sigma, band = make_multiband_lc(band_labels, rng=rng) print(f"Total observations: {len(t)}") print(f"Band counts: { {b: int(np.sum(band == b)) for b in band_labels} }")

Total observations: 240
Band counts: {'g': 80, 'r': 80, 'i': 80}

2. Per-band feature extraction¶

Pass bands= to any feature constructor to enter multiband mode. The output array has len(bands) * n_features_per_band elements.

Band labels can be strings or integer IDs (any NumPy integer dtype). Integer IDs are useful when survey data stores passbands as filter numbers (e.g. ZTF: 1 = g, 2 = r, 3 = i) and are slightly faster because they avoid string hashing during band dispatch.

In [4]:

# Amplitude per band — 3 values (one per passband)
noband_amp = licu.Amplitude()
noband_result = noband_amp(t, m, sigma)

multiband_amp = licu.Amplitude(bands=band_labels)
multiband_result = multiband_amp(t, m, sigma, band)

print(f"All-band amplitude: {noband_result[0]:.3f} mag")
for name, value in zip(multiband_amp.names, multiband_result):
    print(f"{name.removeprefix("amplitude_")}-band amplitude: {value:.3f} mag")

# Amplitude per band — 3 values (one per passband) noband_amp = licu.Amplitude() noband_result = noband_amp(t, m, sigma) multiband_amp = licu.Amplitude(bands=band_labels) multiband_result = multiband_amp(t, m, sigma, band) print(f"All-band amplitude: {noband_result[0]:.3f} mag") for name, value in zip(multiband_amp.names, multiband_result): print(f"{name.removeprefix("amplitude_")}-band amplitude: {value:.3f} mag")

All-band amplitude: 0.822 mag
g-band amplitude: 0.757 mag
r-band amplitude: 0.821 mag
i-band amplitude: 0.662 mag

In [5]:

# LinearFit per band using integer band IDs — no code changes needed in make_multiband_lc
# numpy infers dtype=int64 from integer labels, which the feature accepts directly
int_band_labels = [1, 2, 3]  # e.g. ZTF filter IDs: 1=g, 2=r, 3=i
t_int, m_int, sigma_int, band_int = make_multiband_lc(int_band_labels, rng=rng)

lf_int = licu.LinearFit(bands=int_band_labels)
result_int = lf_int(t_int, m_int, sigma_int, band_int, sorted=True)

print("Feature names with integer band IDs:")
for name, val in zip(lf_int.names, result_int):
    print(f"  {name:35s} = {val:.4f}")

# LinearFit per band using integer band IDs — no code changes needed in make_multiband_lc # numpy infers dtype=int64 from integer labels, which the feature accepts directly int_band_labels = [1, 2, 3] # e.g. ZTF filter IDs: 1=g, 2=r, 3=i t_int, m_int, sigma_int, band_int = make_multiband_lc(int_band_labels, rng=rng) lf_int = licu.LinearFit(bands=int_band_labels) result_int = lf_int(t_int, m_int, sigma_int, band_int, sorted=True) print("Feature names with integer band IDs:") for name, val in zip(lf_int.names, result_int): print(f" {name:35s} = {val:.4f}")

Feature names with integer band IDs:
  linear_fit_slope_1                  = 0.0011
  linear_fit_slope_sigma_1            = 0.0002
  linear_fit_reduced_chi2_1           = 35.8685
  linear_fit_slope_2                  = 0.0014
  linear_fit_slope_sigma_2            = 0.0002
  linear_fit_reduced_chi2_2           = 37.2112
  linear_fit_slope_3                  = 0.0009
  linear_fit_slope_sigma_3            = 0.0002
  linear_fit_reduced_chi2_3           = 32.5998

3. Pure-multiband color features¶

Some features are inherently multiband — they always require bands and have no single-band mode. These compute cross-band statistics rather than per-band statistics.

Feature	Constructor	Output
ColorOfMedian(bands)	exactly 2 bands	median(band[0]) − median(band[1])
ColorOfMaximum(bands)	exactly 2 bands	max(band[0]) − max(band[1])
ColorOfMinimum(bands)	exactly 2 bands	min(band[0]) − min(band[1])
ColorSpread(bands)	≥2 bands	population std dev of per-band weighted mean magnitudes

In [6]:

# g − r median color: filter to configured bands or it raises on unknown passbands
gr = band != "i"
com = licu.ColorOfMedian(["g", "r"])
print(com.names[0], "=", com(t[gr], m[gr], sigma[gr], band[gr], sorted=True)[0].round(4), "mag")

# spread of weighted-mean magnitudes across all three bands
cs = licu.ColorSpread(["g", "r", "i"])
print(cs.names[0], "=", cs(t, m, sigma, band, sorted=True)[0].round(4), "mag")

# g − r median color: filter to configured bands or it raises on unknown passbands gr = band != "i" com = licu.ColorOfMedian(["g", "r"]) print(com.names[0], "=", com(t[gr], m[gr], sigma[gr], band[gr], sorted=True)[0].round(4), "mag") # spread of weighted-mean magnitudes across all three bands cs = licu.ColorSpread(["g", "r", "i"]) print(cs.names[0], "=", cs(t, m, sigma, band, sorted=True)[0].round(4), "mag")

color_median_g_r = -0.0454 mag
color_spread = 0.0181 mag

4. Mixed single-band + multiband Extractor¶

Extractor accepts any combination of single-band and multiband features. Single-band features receive the full light curve; multiband features receive per-band splits. Output names and values are concatenated in declaration order.

ColorOfMedian (a pure-multiband feature from section 3) fits naturally here alongside per-band and whole-curve features:

In [7]:

ext = licu.Extractor(
    licu.Amplitude(bands=band_labels),  # 3 values — per band
    licu.WeightedMean(bands=band_labels),  # 3 values — per band
    licu.ColorOfMedian(band_labels[:2]),  # 1 value  — g − r median color
    licu.ReducedChi2(),  # 1 value  — whole light curve (single-band)
)

result_ext = ext(t, m, sigma, band=band, sorted=True)

print("Feature names and values:")
for name, val in zip(ext.names, result_ext):
    print(f"  {name:35s} = {val:.4f}")

ext = licu.Extractor( licu.Amplitude(bands=band_labels), # 3 values — per band licu.WeightedMean(bands=band_labels), # 3 values — per band licu.ColorOfMedian(band_labels[:2]), # 1 value — g − r median color licu.ReducedChi2(), # 1 value — whole light curve (single-band) ) result_ext = ext(t, m, sigma, band=band, sorted=True) print("Feature names and values:") for name, val in zip(ext.names, result_ext): print(f" {name:35s} = {val:.4f}")

Feature names and values:
  amplitude_g                         = 0.7569
  amplitude_r                         = 0.8208
  amplitude_i                         = 0.6616
  weighted_mean_g                     = 14.9630
  weighted_mean_r                     = 15.0072
  weighted_mean_i                     = 14.9868
  color_median_g_r                    = -0.0454
  chi2                                = 35.1402

5a. Multiband Bins¶

Bins bins observations by time window before evaluating inner features. Add bands= to apply the same binning independently per passband.

In [8]:

bins_mb = licu.Bins(
    [licu.Amplitude(), licu.Mean()],
    window=5.0, offset=0.0,
    bands=["g", "r"],
)

t2, m2, s2, b2 = make_multiband_lc(["g", "r"], n_per_band=100, rng=rng)
result_bins = bins_mb(t2, m2, s2, band=b2)

print("Bins multiband names:", bins_mb.names)
print("Values:              ", result_bins)

bins_mb = licu.Bins( [licu.Amplitude(), licu.Mean()], window=5.0, offset=0.0, bands=["g", "r"], ) t2, m2, s2, b2 = make_multiband_lc(["g", "r"], n_per_band=100, rng=rng) result_bins = bins_mb(t2, m2, s2, band=b2) print("Bins multiband names:", bins_mb.names) print("Values: ", result_bins)

Bins multiband names: ['bins_window5.0_offset0.0_amplitude_g', 'bins_window5.0_offset0.0_amplitude_r', 'bins_window5.0_offset0.0_mean_g', 'bins_window5.0_offset0.0_mean_r']
Values:               [ 0.26937461  0.25340544 15.02650096 15.00719902]

5b. Multiband Periodogram¶

Periodogram supports multiband mode via MultiColorPeriodogram, which finds a single best period that fits all passbands simultaneously. Use multiband_normalization='chi2' for the chi-squared-weighted variant.

In [9]:

pg = licu.Periodogram(peaks=2, bands=["g", "r"], multiband_normalization="chi2")

t_pg, m_pg, s_pg, b_pg = make_multiband_lc(["g", "r"], n_per_band=150, rng=rng)
result_pg = pg(t_pg, m_pg, s_pg, band=b_pg)

print("Periodogram names:", pg.names)
print("Best period (days):", 1.0 / result_pg[0] if result_pg[0] > 0 else "N/A")

pg = licu.Periodogram(peaks=2, bands=["g", "r"], multiband_normalization="chi2") t_pg, m_pg, s_pg, b_pg = make_multiband_lc(["g", "r"], n_per_band=150, rng=rng) result_pg = pg(t_pg, m_pg, s_pg, band=b_pg) print("Periodogram names:", pg.names) print("Best period (days):", 1.0 / result_pg[0] if result_pg[0] > 0 else "N/A")

Periodogram names: ['multicolor_periodogram_period_0', 'multicolor_periodogram_period_s_to_n_0', 'multicolor_periodogram_period_1', 'multicolor_periodogram_period_s_to_n_1']
Best period (days): 0.6378936694123972

6. Batch processing with .many()¶

.many() processes a list of light curves in parallel. For multiband features each element must be a four-tuple (t, m, sigma, band).

In [10]:

feature = licu.Amplitude(bands=["g", "r"])

# Generate 500 two-band light curves
n_lcs = 500
lcs = [make_multiband_lc(["g", "r"], n_per_band=60, rng=rng) for _ in range(n_lcs)]

# .many() returns shape (n_lcs, n_features)
results = feature.many(lcs, sorted=True)
print(f"Shape: {results.shape}")
print(f"Mean amplitude_g: {results[:, 0].mean():.4f}")
print(f"Mean amplitude_r: {results[:, 1].mean():.4f}")

feature = licu.Amplitude(bands=["g", "r"]) # Generate 500 two-band light curves n_lcs = 500 lcs = [make_multiband_lc(["g", "r"], n_per_band=60, rng=rng) for _ in range(n_lcs)] # .many() returns shape (n_lcs, n_features) results = feature.many(lcs, sorted=True) print(f"Shape: {results.shape}") print(f"Mean amplitude_g: {results[:, 0].mean():.4f}") print(f"Mean amplitude_r: {results[:, 1].mean():.4f}")

Shape: (500, 2)
Mean amplitude_g: 0.7023
Mean amplitude_r: 0.7087

7. Arrow input for multiband .many()¶

For zero-copy data exchange from polars / pyarrow / nanoarrow, pass an Arrow List<Struct<...>> array and include "band" in the arrow_fields dict.

See the batch processing tutorial for more complete examples including nested-pandas and Polars.

In [11]:

import numpy.testing as npt
import pyarrow as pa


def make_arrow_lcs(lcs):
    """Convert list of (t, m, sigma, band) tuples to a pyarrow List<Struct>."""
    struct_type = pa.struct([
        ("t", pa.float64()),
        ("m", pa.float64()),
        ("sigma", pa.float64()),
        ("band", pa.utf8()),
    ])
    rows_per_lc = []
    for t_i, m_i, s_i, b_i in lcs:
        rows_per_lc.append([
            {"t": float(t_i[j]), "m": float(m_i[j]),
             "sigma": float(s_i[j]), "band": str(b_i[j])}
            for j in range(len(t_i))
        ])
    return pa.array(rows_per_lc, type=pa.list_(struct_type))


arrow_arr = make_arrow_lcs(lcs[:10])

result_arrow = feature.many(
    arrow_arr,
    sorted=True,
    arrow_fields={"t": "t", "m": "m", "sigma": "sigma", "band": "band"},
)
print(f"Arrow result shape: {result_arrow.shape}")

# Verify against list input
result_list = feature.many(lcs[:10], sorted=True)

npt.assert_array_equal(result_list, result_arrow)
print("Arrow and list results match ✓")

import numpy.testing as npt import pyarrow as pa def make_arrow_lcs(lcs): """Convert list of (t, m, sigma, band) tuples to a pyarrow List<Struct>.""" struct_type = pa.struct([ ("t", pa.float64()), ("m", pa.float64()), ("sigma", pa.float64()), ("band", pa.utf8()), ]) rows_per_lc = [] for t_i, m_i, s_i, b_i in lcs: rows_per_lc.append([ {"t": float(t_i[j]), "m": float(m_i[j]), "sigma": float(s_i[j]), "band": str(b_i[j])} for j in range(len(t_i)) ]) return pa.array(rows_per_lc, type=pa.list_(struct_type)) arrow_arr = make_arrow_lcs(lcs[:10]) result_arrow = feature.many( arrow_arr, sorted=True, arrow_fields={"t": "t", "m": "m", "sigma": "sigma", "band": "band"}, ) print(f"Arrow result shape: {result_arrow.shape}") # Verify against list input result_list = feature.many(lcs[:10], sorted=True) npt.assert_array_equal(result_list, result_arrow) print("Arrow and list results match ✓")

Arrow result shape: (10, 2)
Arrow and list results match ✓

Summary¶

Feature	Constructor	Output
Any Rust feature	Feature(bands=[...])	n_bands × k values per call
Any Rust feature (integer IDs)	Feature(bands=[0, 1, 2])	same, slightly faster — band arrays must be integer dtype
Extractor	Mix freely	Single-band and multiband combined
ColorOfMedian / ColorOfMaximum / ColorOfMinimum	cls(["g", "r"]) or cls([0, 1])	1 value — cross-band statistic
ColorSpread	ColorSpread(["g", "r", "i"]) or ColorSpread([0, 1, 2])	1 value — std dev of band means
Bins	Bins([...], window=…, bands=[...])	Per-band binned features
Periodogram	Periodogram(bands=[...])	Joint period, n_peaks values
.many()	feature.many([(t, m, sigma, band), ...])	(n_lcs, n_features) array
Arrow .many()	add arrow_fields={..., "band": "..."}	Same, zero-copy