camdl Data Model

How external data enters a .camdl model, how CSV columns map to DSL objects, and how dimensions are derived from data.


Long-format data

All external data in camdl is long-format: one row per observation, columns are either index columns (keys into dimensions) or value columns (the numbers the model uses).

Each row is one fact. Each column has one role. This is the tidy data principle applied to compartmental models.

Example: population by district

patch   population
kano_dala   485000
borno_maiduguri 345000
borno_gwoza 78000

One index column (patch), one value column (population). Each row says: “the population of this patch is this number.”

Example: spatial adjacency

src dst weight
kano_dala   borno_maiduguri 0.10
kano_dala   kano_fagge  0.15
borno_maiduguri kano_dala   0.10

Two index columns (src, dst), one value column (weight). Pairs not listed have weight 0 (sparse — most patches aren’t adjacent).

Example: SIA campaign records

patch   day coverage
kano_dala   180 0.82
kano_dala   365 0.75
kano_dala   540 0.71
borno_maiduguri 182 0.68
borno_maiduguri 370 0.55
borno_gwoza 190 0.45

Two index columns (patch, day), one value column (coverage). Kano has 3 rows (3 campaigns), Gwoza has 1. Patches with no campaigns have no rows — coverage is 0 by default (no event).

The pattern: every table is (index₁, index₂, ...) → value. The type signature in the DSL declares which columns are indices (by naming dimensions) and which are values (everything remaining).


The type signature is the mapping

A table declaration maps CSV structure to DSL objects:

name : dim₁ × dim₂ × ... = read("file.tsv")

The dimensions map positionally to index columns in the CSV. The remaining column(s) are the value(s). A CSV with n_dims + 1 columns has exactly one value column (the last).

tables {
  pop : patch = read("data/pop.tsv")
}

The type says 1 dimension (patch). The CSV has 2 columns. Positional mapping:

Column 1 (patch)      → patch dimension index
Column 2 (population) → table value

kano_dala   485000    → pop[kano_dala] = 485000
borno_gwoza  78000    → pop[borno_gwoza] = 78000

For a 2D table:

tables {
  adj : patch × patch = read("data/adj.tsv", default = 0.0)
}

The type says 2 dimensions. The CSV has 3 columns:

Column 1 (src)    → first patch dimension
Column 2 (dst)    → second patch dimension
Column 3 (weight) → table value

kano_dala   borno_maiduguri   0.10  → adj[kano_dala, borno_maiduguri] = 0.10

Column names in the CSV are for human readability. The compiler uses positional mapping from the type signature only. (It does require a header row and skips it.)


dimensions {}: where dimension levels come from

A dimension is a named set of levels — the valid values for that index. age has levels {under5, over5}. patch has levels {kano_dala, borno_maiduguri, ...}.

Levels come from the dimensions {} block, which is declared at the top level of the file. Two forms are supported:

Inline — small structural dimensions:

dimensions {
  age = [under5, over5]
}

Data-derived — large or data-driven dimensions read from a file column:

dimensions {
  patch = read("data/pop.tsv", column = "patch")
}

The read(file, column = "col") form reads the named column, collects unique values in first-occurrence order, and those become the levels of the dimension. All tables referencing patch validate against these derived levels.

After declaring dimensions, stratify(by = X) applies them to compartments:

dimensions {
  age   = [under5, over5]
  patch = read("data/pop.tsv", column = "patch")
}

stratify(by = age)
stratify(by = patch)

tables {
  pop : patch = read("data/pop.tsv")
  adj : patch × patch = read("data/adj.tsv", default = 0.0)
}

Rules

1. Each dimension is defined exactly once.

# OK: patch defined once, referenced by two tables
dimensions { patch = read("data/pop.tsv", column = "patch") }
tables {
  pop : patch = read("data/pop.tsv")
  adj : patch × patch = read("data/adj.tsv", default = 0.0)
}

# ERROR: "dimension 'patch' declared twice in dimensions block"
dimensions {
  patch = read("data/pop.tsv", column = "patch")
  patch = read("data/pop2.tsv", column = "patch")
}

2. Validation is exhaustive and strict.

When a table references a known dimension, the compiler checks the CSV column values against the dimension’s known levels:

  • CSV value not in dimension → error: "'kano_dala_north' in column 1 of adj.tsv is not a valid 'patch' level"
  • Dimension level not in CSV → depends on default:
    • With default = 0.0 → missing cells filled with default (sparse)
    • Without default → error: "patch level 'borno_abadam' has no entry in adj.tsv" (dense — every level must appear)

3. Dimension levels are exhaustively validated.

If stratify(by = patch) appears, the levels from the dimensions {} block become the compartment expansion set. The compiler ensures every level produces a valid identifier (no spaces, no special characters beyond _).

4. Bare dimension names must already be declared.

Unknown bare name in a type signature → error: "unknown dimension 'ptach' — did you mean 'patch'?" (with Levenshtein suggestion). This catches typos that implicit derivation would silently accept.

5. Inline [...] and read(...) are mutually exclusive for the same dim.

# OK: inline levels
dimensions { age = [under5, over5] }

# OK: data-derived levels
dimensions { patch = read("data/pop.tsv", column = "patch") }

# ERROR: "dimension 'age' declared with both inline levels and read()"
dimensions { age = [under5, over5] }
# (cannot also have age = read(...) in the same block)

Two ways to specify dimensions

Source Declaration Stratifies compartments? Levels from
Inline dimensions { age = [under5, over5] } + stratify(by = age) yes inline list
Data, stratified dimensions { patch = read(f, column="c") } + stratify(by = patch) yes CSV column
Data, index-only dimensions { sia_time = read(f, column="c") } (no stratify) no CSV column

All three produce the same thing internally: a named set of levels. They differ in: (a) whether the dimension stratifies compartments, and (b) where the levels come from.

stratify(by = X) declares the role: “this dimension applies to compartments.” The dimensions {} block provides the membership. Role and membership are orthogonal — you can declare a dimension from data without stratifying (e.g., sia_time = read(...) with no stratify), or stratify with inline levels.


Multiple value columns

When a CSV has more columns than n_dims + 1, there are multiple value columns. List multiple table names on the left of ::

tables {
  pop, sex_ratio : patch = read("data/demographics.tsv")
}
patch   population  sex_ratio
kano_dala   485000  0.51
borno_maiduguri 345000  0.49

One index column, two value columns. Creates two tables with the same index: pop[kano_dala] = 485000, sex_ratio[kano_dala] = 0.51.

Value columns map positionally to the names on the left of :. The compiler validates: number of names must match number of non-index columns.

All values are f64. Table values are always numeric. Dimension levels are strings (or numeric-coercible strings). These are different things: levels are the index space, values are what the table stores.


Sparsity and defaults

default = 0.0 fills index combinations with no CSV row:

dimensions {
  sia_time = read("data/sia.tsv", column = "day")
}

tables {
  adj     : patch × patch          = read("data/adj.tsv", default = 0.0)
  sia_cov : patch × sia_time       = read("data/sia.tsv", default = 0.0)
}
Table What default means
adj, default = 0.0 Non-adjacent patches have zero coupling
sia_cov, default = 0.0 No campaign event → zero coverage

Without default, the table is dense: every index combination MUST have a row. Missing rows → compile error. Dense is the right default for things like population — every patch must have a value.

Future: default = missing. For observation data with genuinely absent measurements. The inference engine conditions on present values and skips missing ones. Same table structure, different default semantics. Not needed for forward simulation.


SIA campaigns: the complete pattern

This is where the data model earns its keep. Campaign data is ragged (different patches get different numbers of rounds), indexed by time (which is continuous in the model), and has per-event values (coverage).

The data:

patch   day coverage
kano_dala   180 0.82
kano_dala   365 0.75
kano_dala   540 0.71
borno_maiduguri 182 0.68
borno_maiduguri 370 0.55
borno_gwoza 190 0.45

The model:

dimensions {
  sia_time = read("data/sia.tsv", column = "day")
}

tables {
  sia_cov : patch × sia_time = read("data/sia.tsv", default = 0.0)
}

interventions {
  sia[p in patch, t in sia_time] : transfer(
    fraction = vacc_eff * sia_cov[p, t],
    from = S[under5, p], to = V[under5, p],
    at = t
  ) where sia_cov[p, t] > 0
}

What the compiler does:

  1. Reads sia.tsv. Column “day” → sia_time levels (unique values: [180, 182, 190, 365, 370, 540]). Column 1 → patch (known, validated). Column 3 → value.

  2. Builds sia_cov as sparse 238 × 6 table. Most cells are 0.0.

  3. Expands sia[p in patch, t in sia_time]: Cartesian product of 238 × 6 = 1,428 candidate interventions.

  4. Applies where sia_cov[p, t] > 0: evaluates at compile time (table values are known). ~1,422 zero-coverage cells eliminated. 6 survive.

  5. Each surviving intervention:

    • sia_kano_dala_180: at = 180.0, fraction = vacc_eff × 0.82
    • sia_kano_dala_365: at = 365.0, fraction = vacc_eff × 0.75
    • etc.
  6. IR contains 6 flat interventions. No dimension metadata, no table references. Concrete names, times, expressions.

Numeric dimension levels. sia_time has levels [180, 182, 190, ...]. These are numeric strings from the CSV. They coerce to f64 in expression contexts like at = t. The compiler tracks which dimensions have all-numeric levels and allows this coercion. Non-numeric levels (at = child) would be a type error.

Two campaign data sources. If routine and outbreak campaigns are separate:

dimensions {
  routine_time  = read("data/routine.tsv",  column = "day")
  outbreak_time = read("data/outbreak.tsv", column = "day")
}

tables {
  routine_cov  : patch × routine_time  = read("data/routine.tsv",  default = 0.0)
  outbreak_cov : patch × outbreak_time = read("data/outbreak.tsv", default = 0.0)
}

interventions {
  sia_routine[p in patch, t in routine_time] : transfer(
    fraction = vacc_eff * routine_cov[p, t],
    from = S[under5, p], to = V[under5, p],
    at = t
  ) where routine_cov[p, t] > 0

  sia_outbreak[p in patch, t in outbreak_time] : transfer(
    fraction = vacc_eff * outbreak_cov[p, t],
    from = S[under5, p], to = V[under5, p],
    at = t
  ) where outbreak_cov[p, t] > 0
}

scenarios {
  routine_only  { enable = [sia_routine] }
  full_response { enable = [sia_routine, sia_outbreak] }
}

Two dimensions, two tables, two intervention families. No union, no intersection, no ambiguity. The scenario system composes them.


Time: time_unit and the interpolation boundary

What time_unit does today: normalizes unit literals at compile time. 5 'years with time_unit = 'days becomes 5 × 365.25 = 1826.25. All rate expressions and schedule times are in the declared time unit after compilation. The runtime operates in a single unit — no conversions at simulation time.

The convention: all numeric time values in CSVs are in the model’s time_unit. Campaign day 180 means day 180. Climate week 26 means… day 26? No — the user converts. If time_unit = 'days and the CSV has weekly data, the time column should be 7, 14, 21, ... (days), not 1, 2, 3, ... (weeks). The DSL doesn’t do unit conversion on numeric CSV columns.

ISO-date time columns, however, are converted. If the time column holds ISO dates (2020-03-15, …) instead of numbers, camdl auto-detects this and converts each to internal time via the model’s origin + time_unit — so dated surveillance data loads directly, no hand-anchoring. Numeric columns stay untouched (the convention above). See docs/dates.md for the full rules (auto-detection, --time-format, timezones, output rendering).

Where this matters for at = t: when sia_time has levels [180, 365, 540], at = t produces fire times 180.0, 365.0, 540.0 in the model’s time unit (days). The runtime schedules interventions at exactly these times. No interpolation needed — discrete events at exact times.

Where interpolation IS needed: continuous time-varying covariates.

Climate data, seasonal contact rates, rainfall — these vary continuously but are measured at discrete intervals. The table stores patch × climate_week → temperature, but the Gillespie algorithm asks “what’s the temperature at t = 14.73 days?”

This is NOT a table lookup problem — it’s a time function problem. Tables are compile-time lookups. Time-varying covariates need runtime interpolation.

Per-stratum data-driven forcing: an indexed interpolated time function.

A forcing indexed by a dimension reads each stratum’s own series from a long-format file — one row per (stratum, time) — and interpolates it at runtime. The file has a stratum-id column (key_col), a time column, and a value column:

# data/temperature.tsv
patch     week   temp
kano      0      27.1
kano      28     29.8
katsina   0      26.4
katsina   28     28.9
forcing {
  temperature[p in patch] : interpolated 'ratio {
    data      = "data/temperature.tsv"
    key_col   = "patch"      # the column identifying the stratum
    time_col  = "week"       # numeric time in the model's time_unit, or an ISO date
    value_col = "temp"
    method    = linear   # bare enum: linear | constant | spline
  }
}

transitions {
  infection[p in patch] : S[p] --> I[p]
    @ beta * (1 + alpha * temperature[p]) * S[p] * I[p] / N[p]
}

This expands to one interpolated forcing per patch — temperature_kano, temperature_katsina, … — each carrying only that patch’s rows. At runtime temperature[kano] interpolates the Kano series at the current time t. The stratum a row belongs to is its key_col value, matched to the forcing’s index level; the binder name (p) and the column name (patch) need not coincide, and a level with no matching rows is a compile error (E227), not a silent zero-forcing.

This is a natural extension, not a new concept: time functions already do runtime interpolation, and indexing one by a dimension is the same [p in patch] pattern used everywhere else.

Two data sources. The form above reads its own long-format file (data = …). Alternatively the series can come from a table already in the central tables {} block — the same road population, contact matrices, and adjacency travel — with table = and time_dim =:

tables {
  temp_data : patch × climate_week = read("data/temperature.tsv")
}

forcing {
  temperature[p in patch] : interpolated 'ratio {
    table    = temp_data     # the source table
    time_dim = climate_week  # which dimension is the time axis (its levels are the knot times)
    method   = linear
  }
}

The forcing’s own index ([p in patch]) binds the stratum; time_dim names the time axis; every table dimension must be either indexed by the forcing or the time_dim. A dimension left unaccounted for is a compile error (E229) that names it — never a silent guess:

# clim : patch × season × week, but the forcing only indexes patch —
temperature[p in patch] : interpolated 'ratio { table = clim  time_dim = week }
# E229: table 'clim' has dimension 'season' that is neither indexed by the
#   forcing nor the time axis — index it ([p in patch, s in season]) or
#   aggregate it out.

Use data = for a forcing whose data is a standalone long-format file; use table = when the data is (or should be) a tables {} matrix — reused elsewhere in the model, or when the forcing is stratified over more than one dimension (a single data = file has one key column, whereas table = composes to [p in patch, s in season] + time_dim = week). The knots are sorted by their numeric time, so the dimension’s declaration order does not matter.

For the immediate Nigeria model: per-patch sinusoidal approximation via parameters works:

forcing {
  seasonal = sinusoidal(amplitude = 1.0, period = 365.25 'days,
    phase = phi_season, baseline = 1.0)
}

transitions {
  infection[p in patch] : S[p] --> I[p]
    @ beta[p] * seasonal * S[p] * I[p] / N[p]
}

Shared waveform, per-patch R0 absorbs amplitude differences. Good enough for the forward generative model.


Stress test: every epi data type

Data type Type signature Works?
Population patch → f64 ✅
Population × age patch × age → f64 ✅
Contact matrix age × age → f64 ✅
Spatial adjacency (sparse) patch × patch → f64 ✅ default = 0.0
SIA campaigns (ragged) patch × sia_time → f64 ✅ where > 0 filter
Demographics (multi-value) patch → (f64, f64) ✅ pop, sex_ratio : patch = read(...)
Seroprevalence patch × age → f64 ✅
Environmental covariates patch → f64 ✅
Routine immunization patch × age × ri_month → f64 ✅ 3D table
Detection probability patch → f64 ✅
Genetic distances strain × strain → f64 ✅
Climate / time-varying spatial patch × week → f64 ✅ indexed forcing (long-format data =)
Case data (for fitting) observation data separate pipeline (observations block)

Per-patch time-varying covariates load as an indexed interpolated forcing (a long-format file, one series per stratum — see above). Everything else works with dims → scalar tables, the dimensions {} block, and the existing [i in dim] iteration.


File format rules

Extension determines separator:

  • .csv → comma-separated
  • .tsv → tab-separated
  • Other → error: "unrecognized extension '.txt'; use .csv or .tsv"

Validation:

  • Compiler checks that the detected separator is actually present
  • .csv file with tabs but no commas → error: "file 'data.csv' has .csv extension but appears tab-separated; rename to .tsv"
  • First row is always a header (required)
  • No-header detection: if first row parses as all-numeric when subsequent rows are also all-numeric, warn: "first row of 'data.tsv' looks like data, not a header; add a header row"
  • Comment lines: off by default. read(..., comment = "#") enables skipping lines starting with #

What a row covers: covers and the window columns

A date in an observation file does not say what it means. 2026-07-08 is the day 8 July in a file labelled by the day the count describes, the week 8–14 July in a week-starting file, and the week 2–8 July in a week-ending one. The same cell, three different spans. No rule gets all three right, so an incidence stream states its period rather than having one inferred from row spacing.

A stream that reads an instant — a compartment reading, a proportion prevalence(of = …, among = …), a serosurvey — needs none of this. It has a time column and nothing else; there is no window to state, and declaring one is an error.

Five ways to say it

Four forms take the stream’s own time column and a rule. With D the label in a row’s time column, and the last column showing D = 11 July, seven days:

you write the row covers resolved for
covers = day(time) [D, D+1 day) [11 Jul, 12 Jul) daily counts, labelled by the day
covers = starting_on(time, 7 'days) [D, D+7 days) [11 Jul, 18 Jul) a week labelled by its first day
covers = ending_on(time, 7 'days) [D−6 days, D+1 day) [5 Jul, 12 Jul) — 11th included “week ending D” — ISO/MMWR
covers = closing_at(time, 7 'days) [D−7 days, D) [4 Jul, 11 Jul) — 11th excluded a label that is the closing boundary

Every window camdl writes or reads is half-open: [start, stop) includes start and excludes stop. Two consequences make the notation worth reading carefully. The width is exactly stop − start, with no off-by-one to remember, and two adjacent rows share a boundary value without overlapping — the day ending one window is the day opening the next, counted once. (Readers coming from genomics will recognise the convention: a BED interval spans [start, end) for the same two reasons, and end − start is likewise the length.1)

ending_on and closing_at are near-synonyms in English and differ by a whole day; that is why they are adjacent rows with the interval written out. “Week ending Saturday 11 July” means 11 July is the last included day, so the span runs [5 Jul, 12 Jul) — it starts six days before the label and ends one day after it. That off-by-one is what ending_on exists to hide. closing_at is for a file whose label is the instant the window closes, with that instant itself in the next row’s window: the convention of pomp’s accumulator variables (the value at t[k+1] is the flow over (t[k], t[k+1])), and of camdl’s own simulate --obs output, whose row at t carries the flow since the previous row. If your file was made by either, write closing_at; if it was typed from a surveillance bulletin labelled by the day or week it reports on, write day, starting_on or ending_on.

The fifth form puts both boundaries in the file, one pair per row:

cases_ituri {
  columns     { onset_from : window_start, onset_stop : window_stop,
                cases_ituri : count }
  projected   = incidence(confirm[ituri])
  cases_ituri ~ neg_binomial(mean = p_report * projected, r = k)
}
onset_from  onset_stop  cases_ituri
2026-07-01  2026-07-08  31
2026-07-08  2026-07-15  44
2026-07-15  2026-07-29  96

The third row covers a fortnight because publication was suspended for a week. That is the case with no expression at all under row-spacing inference — the only way to widen a window was to delete the intervening row, which this document separately forbids.

One set of temporal columns, either kind

columns {} requires exactly one set of temporal columns. That is either a time column or a window_start/window_stop pair — never both, and never half a pair. A stream declaring window columns does not also write covers =; the columns are the declaration.

The stop is the stream’s fit time source. Everything downstream that asks “when is this observation” — output rows, --score-from, a forecast grid — reads the closing boundary, so a windowed stream sits on the axis exactly where an unwindowed one does.

That is a position on the time axis and nothing more. It does not say the observation was reported on that date. camdl has no report-date concept: a column carries a time, a dimension, a value, or a window boundary, and nothing else. Reporting delay is a modelled quantity — a delay or convolution term written into the model and estimated like any other — never something read out of a file’s temporal columns. A row whose window is [8 Jul, 9 Jul) is the count for 8 July whatever date sits in the stop column.

A width that varies row to row

Some files carry their own width — a days_covered column saying how many days each row actually represents, because publication slipped. That is what the window columns are for, and the conversion is one subtraction done once where the file is built. The stream names the two boundary columns and nothing else:

cases {
  columns   { window_start : window_start, window_stop : window_stop,
              cases : count }
  projected = incidence(infection)
  cases     ~ poisson(rate = rho * projected)
}
window_start    window_stop cases
2026-07-07  2026-07-08  12
2026-07-08  2026-07-09  9
2026-07-09  2026-07-14  61

The last row is a five-day total. Note there is deliberately no form that takes the width as a duration column (ending_on(time, days_covered)). It would say exactly what the two columns above say, and to do so it would need its own convention about whether the labelled day is included in the count — which is the off-by-one this whole section exists to remove. One way to state a window, stated explicitly.

A gap between rows

The window columns are also the only form that can leave a span uncovered. If the row 2026-07-14 → 2026-07-15 is simply absent from the file above, the flow over that day belongs to no row and is discarded — nothing is scored against it, and the next row’s bin opens at its own window_start. That is the same reading an NA row gets (see the next section): no likelihood term, and the window is kept.

Under the four uniform forms a gap is refused. Every row’s window follows from its label, so consecutive rows whose windows do not touch can only mean a row is missing — and a missing row is exactly what used to widen the next bin silently. The error names both windows and the uncovered span between them. Keep every scheduled row and write the unobserved one as NA; if one row genuinely covers the whole span, state that with the window columns.

Written out, with 9 July simply absent:

cases {
  columns   { time : time, cases : count }
  covers    = day(time)
  projected = incidence(infection)
  cases     ~ poisson(rate = rho * projected)
}
time    cases
2026-07-07  12
2026-07-08  9
2026-07-10  14

day(time) gives those three rows the windows [7 Jul, 8 Jul), [8 Jul, 9 Jul) and [10 Jul, 11 Jul). Nothing covers 9 July, so the load stops there instead of quietly scoring the last row against two days of flow.

What this changes for an existing file

Nothing about the file; everything about how it is read. A stream that declared nothing was scored over (previous row, this row] — a window nobody wrote down. For a file whose label is the closing boundary — anything simulate --obs wrote, anything derived from pomp — that window was right all along, and covers = closing_at(time, <spacing>) states it; nothing in the fitted numbers moves. For a daily file labelled by the day the count describes, that window read the count one day early; declaring covers = day(time) moves the scoring to the day the label names, which is a real change in the fitted numbers and is the correction, not a regression. The declaration is where you say which kind of file you have; there is no default, because either guess is silently wrong for the other kind.

The declaration also settles where scoring begins. The first period opens where it says; the warm-up before it is simulated but not scored, and there is no separate setting for that boundary — the declaration alone places the first bin. A first period that opens before simulate.from is an error, since that time is never simulated.

A fraction of the window’s events

A column that reports “of this window’s events, how many were of one kind” — community deaths among the week’s classified deaths, positives among tests — is a ratio of two counts accumulated over the same window, and the stream that scores it accumulates both: projected divides one flow sum by another, and the row’s period is declared exactly as for a count. The denominator is a column the file supplies; the model’s own denominator enters through projected, not through n.

comm_frac {
  columns     { week_ending : time, comm_deaths : count, n_deaths : count }
  covers      = ending_on(week_ending, 7 'days)
  projected   = incidence(die_comm) / (incidence(die_comm) + incidence(die_fac))
  comm_deaths ~ binomial(n = n_deaths, p = projected)
}
week_ending comm_deaths n_deaths
2026-07-07  3   5
2026-07-14  9   17
2026-07-21  0   0
2026-07-28  14  41

The row for the week ending 21 July classified no deaths. It is not a hole: n_deaths = 0 is an observation, and it scores exactly 0 whatever fraction the model produces, because a binomial with no trials has one possible outcome. A row with n_deaths > 0 in a window where a trajectory produced no deaths at all is refused for that trajectory — the data recorded classified deaths where the model made none — and the refusal names the argument. That is the model disagreeing with the row, not a data error, so it is not softened.

Observed values must be in the family’s support

Every value in an observation file is read in the support of the likelihood that scores it, and a value outside it is refused when the file is loaded. The refusal names the file, the line, the column and the value; nothing is scored.

column must be
the scored column of poisson, neg_binomial, zero_inflated, binomial, beta_binomial a non-negative integer
the trial count of binomial / beta_binomial, when n names a column (n = tested) a non-negative integer
the scored column of bernoulli 0 or 1
the scored column of normal or beta, and any other column the likelihood reads a finite number

A row that breaks the rule is a data error, so it cannot be scored as something nearby: a -4 read as a count would have to become some count, and any choice (0, 4) is an observation nobody made. Correct the row, or write NA if the value was not observed — a hole is never checked, because it has no value. An integer written with a decimal point (4.0) is the integer.

This is about the file, not the fit. A trajectory that cannot produce a row’s count — a nonzero count in a window where the model made no flow — is the model disagreeing with valid data, and is scored as impossible for that trajectory.

A binomial count cannot exceed its trial count

For binomial and beta_binomial, the scored count is the number of successes out of n trials, so each observed row must have count <= n, and n must be a non-negative whole number. Whether that is checked when the file is loaded depends on what n reads.

n is checked
a column (n = tested), a literal (n = 50), or arithmetic on columns and literals at load, for every row with both present
anything that reads a compartment, a parameter, time, a table or a named constant at run time, per trajectory

At load, a row with count > n, a negative n, or an n that is not a whole number (n = tested / 2 with tested = 91) is refused. The refusal gives the stream, the row, the time and both values; nothing is scored. A row where n is 0 and the count is 0 is valid and scores exactly zero; a positive count against n = 0 is refused. A row with NA in the count or in a column n reads is a hole and is not checked. Wrap a fractional expression in floor() or ceil() to choose the rounding yourself; the likelihood would otherwise round it, and the trial count scored would not be the one written.

At run time, n takes a different value on each trajectory, so count > n is the model disagreeing with valid data, and that trajectory is scored as impossible (log-likelihood -inf), exactly as for any other family. When every particle is impossible at an observation, PGAS initialization and pfilter both say so and name the cause: the observation’s time, the stream, the observed value, the range of what the particles projected, and how many particles each guard rejected, with the count exceeding the trial count as one such guard. pfilter then exits with status 3 (the filter rejected these parameters) and stores no result, instead of printing -inf as the answer; with --replicates, it does so when every replicate lost all its particles, and otherwise reports the mean as -inf with a warning naming the first collapse. A collapse that persists over several consecutive observations stops the filter earlier, with its effective-sample-size collapse error (also status 3), which does not name the guard. An n that evaluates to a negative number or to NaN is not read as zero trials: it is outside the family’s domain and also scores -inf, reported as the n argument being out of domain. A fractional n at run time (an n built from a real-valued compartment) is rounded to the nearest whole number.

camdl survey reads the same n = <column> as camdl fit and applies the same rules. Survey scores every row, so an NA in a column the likelihood reads is refused rather than made a hole.

Missing observations: NA is a hole, not a zero

In an observation file, the token NA in a value column marks a hole: a time at which the stream was scheduled but no value was observed.

A hole is not an observed zero, and the difference is load-bearing:

  • A hole contributes no likelihood term. It is marginalised out — the filter neither rewards nor penalises the model for what it did not see. An observed 0, by contrast, is a real datum saying “we looked and counted none,” and the model is scored against it.
  • A hole still closes an incidence bin. For a stream whose projected is an incidence(...), the accumulator resets on the observation SCHEDULE, not on value presence (the same convention as pomp’s accumvars). So a row with NA ends the interval at its own time, and the NEXT row is scored over the span since that hole — not since the last row carrying a number.

That second property is the reason to write NA rather than delete a row. If a row is dropped from the file entirely, the following row’s exposure window silently widens to span the gap, and its count is then compared against more modelled flow than it represents. Keeping the row as NA preserves the window; deleting it changes the meaning of the row after it.

Both forms are therefore legitimate and mean different things:

you write the model sees
2026-07-14 NA scheduled, unobserved — no term, interval closes
2026-07-14 0 observed zero — scored as a zero count
(no row at 2026-07-14) not scheduled — the next row’s interval spans it

A daily file with one unobserved day, and the stream that reads it:

cases {
  columns   { time : time, cases : count }
  covers    = day(time)
  projected = incidence(infection)
  cases     ~ poisson(rate = rho * projected)
}
time    cases
2026-07-07  12
2026-07-08  NA
2026-07-09  14

The middle row states the window [8 Jul, 9 Jul) and scores nothing in it. The row after it opens at 9 July, where its own declaration puts it, rather than reaching back over the day nobody observed.

Under a declared covers, a hole needs no separate rule. The two properties above were a consequence of inference: the accumulator reset on row spacing, so a hole had to close a bin to stop the next row’s window swallowing it. With the period stated, an NA row is simply a row whose period scores nothing — it opens and closes exactly where it says, and the row after it opens where its own declaration says. The behaviour is unchanged; it stops being a convention.

This also settles a case that had no defined answer before: a row whose window reaches back across an earlier NA row. Under inference that was ambiguous — the hole had already closed the accumulator, so the later row would be scored against less flow than its window claimed, silently. Under declaration it is simply two overlapping periods, which is an error naming both rows. A stream whose rows genuinely overlap is not expressible, because the flow in the shared span would be scored twice.

A warning for downstream consumers. NA is not a null token in every tool. In polars it does not parse as null by default, so one NA turns a numeric value column into a string column; a plotting library then treats those strings as categorical and can draw them at the axis baseline — holes rendering as observed zeros in a figure that looks entirely plausible. Read observation files with NA declared as the null value (in polars, null_values="NA").

Building the file: emit the schedule, not the usable subset

The rule above is stated from the reader’s side. From the writer’s side it inverts into one instruction, and it is the instruction most often got wrong:

Emit a row for every scheduled time. Make the unusable ones NA. Never filter them out.

A pipeline that drops bad days is doing the natural thing — a filter is how you remove data you cannot trust. But for an incidence(...) stream a filter does not remove a day; it merges that day into the next observation’s window, while the surviving row’s count still covers only the one day it was reported for. The window widens; the count does not. (For which calendar day a dated row covers, see camdl docs dates, “What a dated row covers”.)

What that costs, from a real case. A daily lab stream whose unusable days were filtered rather than nulled left the retained rows up to 13 days apart. Each was scored against 13 days of accumulated modelled flow and one day of specimens — a 13x inflation of projected/tests, which refused 23 of 24 chains at initialisation. That refusal is correct: a modelled flow above the specimen count is impossible, not merely a poor fit. But nothing in the diagnostic mentions the filter, and the pipeline that caused it looks like ordinary data cleaning.

Re-indexing onto the full span is the fix, and it is three lines:

span = pl.date_range(out["time"].min(), out["time"].max(),
                     interval="1d", eager=True).alias("time")
full = pl.DataFrame({"time": span}).join(out, on="time", how="left")
full.write_csv(path, separator="\t", null_value="NA")

Every scored bin is then exactly one cadence wide, and the days you could not use contribute no observation at all.

Proportion streams need a denominator on the unscored rows too. For a stream scored as k / n, fill n with 1 — or any positive value — where the observation is NA. The value is never read, because the observation is a hole; the fill exists so a 0 denominator cannot reach a division. Leaving it null or zero makes the file’s validity depend on how the projection happens to be written.

The count you should check. After building, assert that the number of rows equals the number of scheduled times, and that the number of non-null rows equals the number of usable observations. Those two numbers differing by exactly the number of holes is the property this section exists to protect; if the row count equals the usable count, the filter is still there.


Complete Nigeria model with data model

time_unit = 'days

compartments { S, E, I, R, V }

dimensions {
  age      = [under5, over5]
  patch    = read("data/demographics.tsv", column = "patch")
  sia_time = read("data/sia.tsv", column = "day")
}

stratify(by = age)
stratify(by = patch)

parameters {
  sigma    : rate        in [0.001, 1.0]
  gamma    : rate        in [0.001, 1.0]
  omega    : rate        in [0.0001, 0.1]
  kappa    : rate        in [0.0, 0.5]
  vacc_eff : probability in [0.0, 1.0]
  I0       : count       in [1, 10000]
  R0[patch] : positive   in [0.5, 15.0]
}

let beta[p in patch] = R0[p] * gamma
let N[a in age, p in patch] = S[a, p] + E[a, p] + I[a, p] + R[a, p] + V[a, p]

tables {
  # patch levels derived from demographics.tsv (238 LGAs)
  pop, init_sus : patch                = read("data/demographics.tsv")

  adj      : patch × patch             = read("data/adj.tsv", default = 0.0)
  age_frac : age                       = [0.18, 0.82]
  sia_cov  : patch × sia_time          = read("data/sia.tsv", default = 0.0)
}

transitions {
  infection[a in age, p in patch] : S[a, p] --> E[a, p]
    @ beta[p] * S[a, p] * I[a, p] / N[a, p]

  importation[a in age, p in patch, q in patch] : S[a, p] --> E[a, p]
    @ kappa * adj[p, q] * S[a, p] * I[a, q] / N[a, q]
    where p != q

  progression[a in age, p in patch] : E[a, p] --> I[a, p]  @ sigma * E[a, p]
  recovery[a in age, p in patch]    : I[a, p] --> R[a, p]  @ gamma * I[a, p]
  waning[a in age, p in patch]      : V[a, p] --> S[a, p]  @ omega * V[a, p]
}

interventions {
  sia[p in patch, t in sia_time] : transfer(
    fraction = vacc_eff * sia_cov[p, t],
    from = S[under5, p], to = V[under5, p],
    at = t
  ) where sia_cov[p, t] > 0
}

init {
  S[a in age, p in patch] = pop[p] * init_sus[p] * age_frac[a]
  I[under5, borno_damboa] = I0
}

simulate {
  from = 0 'days
  to   = 730 'days
}

scenarios {
  baseline {
    label = "no SIA"
  }
  with_sia {
    label = "with SIA"
    enable = [sia]
  }
}

Summary

Concept Mechanism Example
Small structural dimension dimensions { X = [...] } age = [under5, over5]
Large data-driven dimension dimensions { X = read(f, column="c") } patch = read("pop.tsv", column="patch")
Index-only dimension dimensions { X = ... } (no stratify) sia_time = read("sia.tsv", column="day")
Dimension validation bare X in type signature adj references known patch
Table loading read("file.tsv") positional column → dimension mapping
Sparse tables default = 0.0 missing index combinations → default
Multiple values a, b : dims = read(...) two tables from one CSV
Numeric level coercion at = t where t iterates numeric dim campaign times as floats
Compile-time filtering where expr > 0 skip zero-coverage events
Time-varying spatial indexed forcing, long-format data = temperature[p in patch]

One loader, one table type, one iteration syntax. All external data is dims → scalar. All iteration is [i in dim]. Dimensions are declared in the dimensions {} block, validated everywhere. No union, no intersection, no implicit derivation. Every mapping from CSV to model is explicit and traceable.

Footnotes

  1. UCSC Genome Browser, “BED format”, https://genome.ucsc.edu/FAQ/FAQformat.html#format1.↩︎