Constrained stratified allocation

Distribute a total sample size across strata defined by a single stratification variable, under a fixed total n, target CV, or budget. When the design uses multiple stratification variables (e.g. region and urbanicity), cross them into a single variable beforehand so that each row of frame represents one unique stratum.

Usage

n_alloc(frame, ...)

# Default S3 method
n_alloc(
  frame,
  ...,
  domains = NULL,
  n = NULL,
  cv = NULL,
  budget = NULL,
  measures = NULL,
  targets = NULL,
  psu = NULL,
  certainty_cutoff = NULL,
  n_psu_per_zone = NULL,
  control = NULL,
  objective = NULL,
  alloc = c("neyman", "optimal", "proportional", "power"),
  unit_cost = NULL,
  alpha = 0.05,
  deff = 1,
  resp_rate = 1,
  eligibility_rate = 1,
  eligible_per_unit = 1,
  df = NULL,
  min_n_stratum = NULL,
  min_n_domain = NULL,
  alloc_q = 0.5,
  fpc = c("unit", "stage", "none"),
  plan = NULL
)

# S3 method for class 'svyplan_prec'
n_alloc(frame, ..., n = NULL, cv = NULL, budget = NULL)

Arguments

frame

For the default method: a stratum-level data frame describing the population you want to sample. Each row represents one stratum, a subgroup of the population defined by a stratification variable such as region, age group, or urbanicity. The values in this frame typically come from a census, a population register, or a previous survey. Any stratum table with the columns below works, including a pool summary carried back from an executed sample once the measure columns are added.

When a design stratifies by several variables at once (e.g. region $\times$ urbanicity), cross them into a single variable before calling n_alloc (e.g. with interaction()) so that each row maps to exactly one population cell.

When measures and targets are supplied, frame instead describes atomic allocation strata for joint constrained allocation. It requires unique, non-missing stratum and positive finite N columns. It may contain unit_cost, max_weight, take_all, and complete domain classification columns named by targets$domain. In that mode sd, var, mean and p belong in measures instead. Each row must be homogeneous for every domain classification used by a target. For a fixed-take two-stage design, also supply whole N_psu and n_per_psu, plus positive cost_psu and cost_ssu. A three-stage design additionally requires whole N_ssu and n_per_ssu, and positive cost_tsu. N_psu is the number of PSUs in the stratum. In three-stage mode, N_ssu is the aggregate number of SSUs across all PSUs in the stratum (not the number per PSU). Stage costs replace unit_cost.

Required columns:

N

Number of units (e.g. households, individuals) in each stratum. These are population counts, not sample sizes. Must be positive and finite.

sd or var (classic modes)

A measure of how spread out the variable of interest is within each stratum, required except for a fixed n under alloc = "proportional" or "power". Provide exactly one:

sd: the stratum standard deviation ($\sqrt{\text{variance}}$), or
var: the stratum variance.

Both must be non-negative and finite. A frame with p and neither column takes its variance from p, as measures does in joint allocation, $N p (1 - p) / (N - 1)$ in an element design and $p (1 - p)$ in a multistage one or under fpc = "none". Without sd, var or p, those two rules allocate a fixed n from the stratum sizes alone, and the precision fields of the result are NA because nothing was assessed.

Optional columns:

stratum: A label identifying each stratum (e.g. "Urban", "Rural"). If omitted, row numbers are used. Must be unique, or unique within each domain when domains is set.
Label columns: Any other character, factor or logical column, such as region and residence, is carried into $detail after stratum. Numeric columns are not, so an integer code such as a region number is carried only once converted with factor(). A label column named like a column of $detail is an error.
mean or p: The stratum population mean or proportion of the variable of interest. Required when solving for cv, and by alloc = "power" when sd or var is given, because the coefficient of variation is defined relative to the mean. Use mean for continuous variables and p (in $[0, 1]$) for binary (yes/no) variables.
alloc_measure: Positive finite stratum measure used only by Bankier allocation (alloc = "power"). Defaults to N * abs(mean) (or N * p). Population CV is derived from sd / abs(mean) and requires a nonzero mean in every stratum. Without sd, var or p the CV term drops, the weights are alloc_measure^alloc_q, and alloc_measure defaults to N.
unit_cost: Per-unit interviewing cost in each stratum (positive, finite). Set higher values for strata that are more expensive to reach. Defaults to 1 everywhere (equal cost).
max_weight: Maximum allowed sampling weight $N_h / n_h$. Caps how under-represented a stratum can be. Use NA for strata without a cap. It bounds the gross weight, the units drawn, not the response-adjusted $N_h/(n_hr_h)$, which differs once resp_rate varies by stratum.
min_n_stratum: Minimum sample size of each stratum, NA for none. See the argument of the same name.
n_step: Whole step of each stratum's sample size, NA or 1 for none, so that n_int is a multiple of it. Use it when the sample of a stratum must divide evenly over fieldwork periods or rotation groups, such as 12 PSUs a year for monthly fieldwork. It counts the units n counts, so on a frame of households with n counting PSUs it is a number of PSUs. Bounds move to whole steps. A floor rises to the next step, the population caps a stratum at its last whole step, and a take_all stratum keeps its N. A two-stage frame and joint allocation refuse it. See the section on steps in Details.
take_all: Logical (or 0/1). If TRUE, every unit in the stratum is included, a census stratum. Useful for small strata whose total population is tiny enough to enumerate.
deff, resp_rate: Per-stratum design effect and expected response rate, for designs whose fieldwork differs across strata. NA in either column falls back to the scalar argument. The deff and resp_rate arguments override these columns when given as vectors. In joint allocation a non-missing measures row overrides them instead. A published DEFT is squared first, and the allocation already carries its own weighting loss. See varcomp().

For svyplan_prec objects: a precision result from prec_alloc().

...

Additional arguments passed to methods. Unused arguments are rejected.

domains

Character vector of column names in frame to treat as domain identifiers, or NULL (default) for no domains. All names must exist in frame. Domains define sub-populations that each contain one or more strata. When cv is the target, precision is enforced within every domain (see Details).

n

Total sample size. Specify exactly one of n, cv, or budget.

cv

Target coefficient of variation (relative standard error). For example, cv = 0.05 means the standard error of the estimated population mean or total should be at most 5 percent of the estimate. Requires mean or p in frame. When domain columns are present, this target is enforced in each domain. Specify exactly one of n, cv, or budget.

budget

Total field budget. In the classic modes, specify exactly one of n, cv, or budget, and budget caps the affordable allocation. In joint constrained allocation it must be accompanied by objective. The allocation then minimizes the priority-weighted relative variance of the objective components subject to the budget, every hard target in targets, and the stratum bounds. With psu the budget holds the field design's cost, certainty visits included. Supplying budget without objective in joint mode is an error, because the two readings are incompatible.

measures

Optional long data frame for joint constrained allocation, with one row per required stratum and indicator name. Each row must contain either p in [0, 1], or mean and exactly one of non-negative sd or var. In element-level designs, p implies the finite-population variance N * p * (1 - p) / (N - 1), using its stratum's N > 1, consistently with prec_prop(), and p * (1 - p) under fpc = "none". Explicit var and sd values are used as supplied. Multistage designs, including PSU-register designs, retain the working proportion variance p * (1 - p). Optional row-specific deff, resp_rate and eligible_per_unit values take precedence, and NA falls back to the stratum's frame column of the same name, then to the scalar argument. A published DEFT is squared first, and the allocation already carries its own weighting loss. See varcomp(). Indicator values (p, mean, sd, var) belong here, and frame columns of those names are refused. In fixed-take multistage mode, icc_psu is required and var_ratio_psu defaults to 1. Three-stage mode also requires icc_ssu and derives var_ratio_ssu as var_ratio_psu * (1 - icc_psu) when it is absent, the value for SSUs of equal size or drawn proportional to size (see design_effect()). Stage parameters may instead be supplied as stratum defaults in frame, in which case non-missing measure-row values take precedence. Only rows selected by a target and its domain are value-validated. Unused rows are ignored after their non-empty, unique keys and frame-stratum membership are checked. Must be supplied with targets.

targets

Optional long data frame with one precision requirement per row. It requires indicator name and exactly one positive cv, moe, or rmoe. Optional domain and level select a frame domain (.overall and level = NA denote the whole population), alpha overrides the scalar value for MOE, and constraint supplies a stable unique identifier. Without a level column, a row naming a domain stands for every level of that domain, in the order the levels first appear in frame, and takes generated constraint identifiers. With a level column, every domain row names its level, so a table that targets every region and one named district writes the regions out in full. Must be supplied with measures. Optional when objective and budget are given, which requests the best design a budget can buy with no hard precision requirement.

psu

Optional PSU register for certainty-aware allocation, one row per PSU, with stratum and N, the count of ultimate units in that PSU. The PSU sizes must sum to frame$N within each stratum, and every PSU must hold at least its stratum's n_per_psu unless flagged certainty. Optional psu_id is carried through and optional certainty forces a PSU into the certainty part whatever its size. It is an alternative to the frame column N_psu, which says only how many PSUs a stratum has, and requires the frame column n_per_psu and measures$icc_psu. A stratum with a single PSU has no between-PSU variance, so its icc_psu is taken as 0 and may be NA. A stratum that draws a single PSU from its remainder leaves no variance estimate there, and n_alloc() warns with class svyplan_single_draw. n_psu_per_zone = 2 draws at least two. Joint constrained allocation only, to targets or to an objective within a budget, and two stages only. See the certainty section under Details.

certainty_cutoff

Inclusion probability at which a PSU of the psu register is taken with certainty, in $(0, 1]$. One value applies to every stratum, and one per stratum, as a vector or as a frame column of the same name, relaxes the rule where it is set. The argument overrides the column. NULL (default) uses the column where present and 1 otherwise, which is the rule without a cutoff. Requires psu. See the certainty section under Details.

n_psu_per_zone

NULL (default), 1 or 2, the PSUs drawn from each zone when each stratum's remainder is cut into zones of about equal size. With 2, every stratum that draws from its remainder has an estimable variance. With 1, the variance is estimated by collapsing zones in the groups $psu$.pair fixes before selection. An optional psu$zone_order column sets the order zones are cut in, such as a geographic order, and size, largest first, is the default. Requires psu. See the certainty section under Details.

control

NULL (default) or a list of settings for the search that plans the whole design of a joint allocation, with or without psu, and of the classic cv and budget modes that search (see n_int under Value). max_nodes, one positive whole number, 5,000 by default, is the number of search nodes after which the search stops with the best design found, the status "limit" and a lower bound. It finishes the node it is expanding, so it can open a few nodes past the limit. A larger value can prove a design the default leaves at the limit, at the cost of time. gap, a number in $[0, 1)$, 0 by default, stops the search once the design is within that relative gap of the bound, with the status "gap". A stop at the node limit within 1e-6 of the bound also reads "gap". A mode that does not search refuses control.

objective

Optional estimates whose priority-weighted relative variance is minimized among the allocations that meet targets and the budget. Either a character vector of indicator names (overall domain, equal priority) or a data frame with name plus optional domain, level, non-negative priority, and component columns, using the same indicator-domain identifiers as targets, including the rows without a level that stand for every level of a domain. At least one priority must be positive. Priorities need not be normalized. Requires budget and must be supplied with measures.

alloc

Allocation rule: "neyman" (default), "optimal", "proportional", or "power".

unit_cost

Optional scalar or length-nrow(frame) vector of per-stratum unit costs, overriding frame$unit_cost. Not used with fixed-take multistage joint allocation, which requires stage costs.

alpha

Significance level, default 0.05. alpha = 2 * pnorm(-2) gives the multiplier of 2 that the DHS and MICS manuals use in place of 1.96.

deff

Design effect multiplier (> 0). A scalar applies to every stratum. A length-nrow(frame) vector gives one per stratum and overrides a deff column in frame. A published DEFT is squared first, and the allocation already carries its own weighting loss. See varcomp().

resp_rate

Expected response rate among eligible units, in (0, 1]. Default 1. A scalar applies to every stratum. A length-nrow(frame) vector gives one per stratum and overrides a resp_rate column in frame. See the nonresponse section of svyplan-package for what the rate counts and what this adjustment does and does not claim.

eligibility_rate

Share of issued units in scope, such as dwellings found occupied, in (0, 1]. Default 1. A scalar applies to every stratum. A length-nrow(frame) vector gives one per stratum and overrides an eligibility_rate column in frame. In joint allocation it is a scalar, a frame column or a measures column. Strata are weighted by their population in scope, N * eligibility_rate, while the bounds on the sample stay in frame units. With a PSU register it is a scalar or a frame column, one rate per stratum. See the nonresponse section of svyplan-package.

eligible_per_unit

Eligible units per sampled unit in joint allocation, a positive number, default 1, such as married women per household. It can exceed 1. A measures column of that name sets it per indicator and stratum, and a frame column per stratum, with the precedence resp_rate follows. N, n, bounds, weights and costs stay in sampled units, while every variance term reads N * eligibility_rate * eligible_per_unit eligible units and resp_rate * eligibility_rate * eligible_per_unit respondents per sampled unit. Classic allocation and a psu register refuse it. See n_prop() for the caveat above 1.

df

Degrees of freedom of the variance estimator the allocation will have, typically sampled PSUs minus strata, and available from design_df(). It switches the quantile used to translate a moe target in and to report moe out. A cv target carries no quantile and is unaffected. NULL (default) applies no adjustment.

min_n_stratum

Optional minimum sample size per stratum. One value applies to every stratum. A vector of length nrow(frame), or a frame column of the same name, gives each stratum its own floor, with NA for none, such as a domain that needs a set number of interviews whatever its share of a power allocation. The argument overrides the column. The minimum is a lower bound and not a base to which a share is added. A stratum whose share already exceeds its minimum keeps that share, and the strata held at their minimum leave less for the others. To give every stratum its minimum and share the rest on top of it, allocate n - sum(min_n_stratum) and add the minimums, as in the examples.

min_n_domain

Optional floor on each domain's total sample, in the units n counts, such as a publication minimum of 400 households per province. One number applies to every domain. A data frame with columns domain (the domains column), level and min_n gives floors that differ, and a level it leaves out has none. Needs one domains column and a classic mode. See the domains section in Details.

alloc_q

Bankier allocation exponent, used only when alloc = "power". Numeric scalar in $[0, 1]$. With the default alloc_measure, alloc_q = 1 gives Neyman allocation. At alloc_q = 0 weights depend on population CVs and the response/design-effect adjustment, independently of alloc_measure. This does not generally equalize domain CVs. Default 0.5.

fpc

Which finite population correction the variance carries: "unit" (default) for the ultimate-unit 1 - n / N, on the whole inflated variance in cluster allocation, "stage" for the stage-by-stage correction, exact under the equal-size, equal-take variance-component model cluster allocation assumes, or "none" for no correction, n_cluster()'s default with-replacement model. "stage" and "none" match n_cluster() with the same fpc. "stage" requires N_psu, and N_ssu as well at three stages, and is refused in element allocation. Only "unit" is available with a psu register. See the correction section in Details.

plan

Optional svyplan() object providing design defaults.

Value

A svyplan_n object with type = "alloc" and a stratum-level allocation table in $detail (also available via as.data.frame()), with columns:

stratum, N, sd, unit_cost: Stratum identifiers and inputs carried over from the frame. The frame's label columns follow stratum, and $params$strata_keys names them when together they identify the strata, NULL otherwise. In cluster mode unit_cost is the derived effective per-element cost cost_psu / n_per_psu + cost_ssu (1 when no stage costs were given), while sd stays the input stratum SD.
n: Allocated sample size (continuous).
n_int: Integer allocation. In n mode an element design keeps the requested total (bounded largest-remainder rounding). A cluster design fields whole PSUs at one whole take per stratum, near its cost-optimal value, so its total can miss n by at most half the largest take. The allocation is re-solved within what each stratum can field at its take, and each rounded PSU count then moves by at most one PSU, jointly across strata, toward n. A take is raised above that value only when the strata at their N_psu limits would otherwise fall more than half a take short, and a take risen past N_psu under "optimal" fields every PSU at its nearest whole value. In cv mode with alloc = "optimal" or "neyman", the design is the cheapest whole allocation meeting the target (the smallest for "neyman"), with domains the cheapest in each domain, and in budget mode with "optimal" and no domains the whole allocation with the smallest variance within the budget, each found by a search reported in $operational$status and $optimization$search (with domains, the search of the domain furthest from its bound, each domain's search in $optimization$domains). A cluster design searches whole PSU counts and, under "optimal" with stage costs and N_psu, each take the frame leaves free, at any whole value up to the average PSU size N / N_psu that the stratum's bounds allow. A fixed take, every take under "neyman", and every take without N_psu stays at the whole take nearest its value. The other cv designs round each stratum up so the integer design meets the target, and the other budget designs add units by variance reduction per unit cost so the integer design stays within budget. For an optimal design across domains, use joint allocation with measures and targets. Always inside the integerized bounds (ceiling(.lower), floor(.upper)), except that an N_psu cap bounds the whole PSU count rather than the element count, since the whole take can round up. An error is raised when no integer allocation satisfies them, and when no whole-unit cluster design comes within half a take of n, above or below. That error gives the nearest total the design can field, the smallest when lower bounds and one PSU per stratum force more, the largest when N and N_psu allow less. With an n_step column every size is a whole number of steps, as the section on steps in Details describes, and $detail carries n_step.
weight: Base weight N / n of the issued sample, before nonresponse or calibration adjustment.
n_eff: Effective sample size n * resp_rate / deff.
.lower, .upper: Bounds applied to the stratum (from min_n_stratum, max_weight, take_all, N_psu, or N).
.binding: Whether the allocation sits on one of its bounds.
.bound_source: Which constraint produced the bound the allocation sits on, and NA when it sits on none. "min_n_domain" marks the strata of a domain held at its floor.
.psu_frac: Cluster mode with N_psu only: the share of the stratum's PSU universe the design selects, n_psu / N_psu.
.se, .moe, .rmoe, .cv: Precision of the stratum's own mean under the allocation, on the same scale as mean (or p). Multiply .se by N for the stratum total. .rmoe is .moe as a fraction of the stratum mean. .cv and .rmoe are NA unless the frame supplies mean or p.
.share: The stratum's share of the design variance of the overall mean. Sums to 1 and identifies where precision is actually being bought. It is NA when the design variance is zero, which a census is, since no stratum then holds a share of nothing.
take_all, mean: Present when take-all strata or stratum means were supplied.
n_per_psu, n_psu: Cluster mode only: the continuous per-stratum take and implied number of PSUs (n / n_per_psu). A free take above the cost-optimal one is a stratum fielding all its N_psu PSUs (see the N_psu section in Details).
n_psu_int, n_per_psu_int: Cluster mode only: the whole-unit field design (whole PSUs and whole takes), chosen so that budget designs stay within budget and cv designs meet the target. Here n_int = n_psu_int * n_per_psu_int.
n_psu, n_psu_int, n_per_psu, n_per_ssu: Joint fixed-take multistage mode: continuous and whole first-stage decisions and the fixed later-stage takes. n_per_ssu is present only for three stages. Public n remains in ultimate-unit units, so n = n_psu * n_per_psu for two stages and n = n_psu * n_per_psu * n_per_ssu for three stages. The same exact identities hold for n_int and n_psu_int.
n_psu_certain, n_psu_rest, threshold: With psu only: the number of certainty PSUs, the number of PSUs left in the remainder, and the certainty threshold certainty_cutoff * n_per_psu / f, at the sampling fraction f whose classification and draw give the planned design.
n_certain_int, n_psu_draw: With psu only: the units the operational design fields in the certainty PSUs, and the number of remainder PSUs it draws. Here n_int = n_certain_int + n_psu_draw * n_per_psu.
n_zone: With n_psu_per_zone only: the zones the remainder is cut into, so n_psu_draw = n_zone * n_psu_per_zone.

The result also carries an $operational list describing the integer field design. In joint mode it contains total n, variable field cost, the full target table in constraints, all_pass, and status, "optimal" when the search proved the design the best, "gap" when it stopped within the gap of control, or within 1e-6 at the node limit, and "limit" when its node limit came first with a larger gap. With psu it adds n_certain and n_psu_draw. In the classic modes it also reports the applicable aggregate se, moe, and cv. Top-level quantities describe the continuous design, while as.integer() returns the operational total.

With psu, the result adds $psu, one row per PSU of the register:

psu_id, stratum, N: Carried over from the register. psu_id is present when the register has it.
certainty: Whether the PSU is in the certainty part.
n_take: The take the operational design fields in the PSU. A certainty PSU carries its own whole take at the stratum rate, capped at its size, and a remainder PSU carries n_per_psu.
.certainty_source: Why a certainty PSU is certainty. "threshold" when its size reaches the threshold, "supplied" when psu$certainty flagged it below the threshold, "operational" when it is below the threshold but the whole-PSU remainder draw of $operational would give it an inclusion probability of one, and "planned" when the planner holds it below the threshold to cut clustering, which it does only under stage costs. NA for a remainder PSU. summary() counts them.
.threshold, .distance: The stratum's threshold, and the relative distance N / .threshold - 1. A PSU with a small .distance is the one whose classification the allocation can tip.
.zone: With n_psu_per_zone only: the remainder PSU's zone within its stratum, NA for a certainty PSU.
.pair: With n_psu_per_zone = 1 only: the variance group of the PSU's zone, numbered across strata, NA outside any zone.

In joint mode, $optimization$search reports how the whole design was chosen:

family: "whole units" without psu, where a stratum's designs are its whole sizes or whole PSU counts. With psu, "voluntary" under stage costs, where a stratum's designs include those holding its largest PSUs in certainty below the threshold, and "closure" under unit cost, where they do not.
status: "optimal", "gap" or "limit", as in $operational.
bound, gap: A lower bound on the cost of any design in the family that meets every target, and the relative gap (cost - bound) / cost, 0 when the design is proved the cheapest. In budget mode the bound is on the weighted relative variance of any design within the budget, and the gap is relative to the design's.
nodes, max_nodes: The search nodes used and the limit, control$max_nodes.
states, kept: The stratum designs listed, and those the bound left in the search.

In fixed-budget objective mode the result adds $objective, one row per component with its priority, relative variance .relvar, equivalent .cv, weighted .contribution, and .share of the total, plus the scalar $objective_value. $operational carries the same two fields for the integer design together with budget_residual. $params$mode is "budget_objective" rather than "targets", and $optimization adds the budget, whether it binds, its residual, the local budget_sensitivity (change in objective value per extra unit of budget), and the epsilon-constraint objective_bound and its multiplier. With psu these describe the continuous relaxation, and $operational the field design the search planned.

Details

Joint constrained allocation

Supplying both measures and targets requests a minimum-cost allocation meeting every indicator-domain precision requirement at once, and adding objective with budget asks instead for the best design that budget can buy. This generalized mode is mutually exclusive with scalar n, cv, budget, domains and an explicit alloc rule, because the target table defines the constraints and the allocation is solved jointly. Its input contract, its multistage variance model, its infeasibility messages and the predict() frontier are on n_alloc-generalized.

For fieldwork, use $detail$n_int, not the analytical continuous n. Fixed-take count columns are numeric vectors whose values are exact whole numbers. The _int suffix denotes operational integer semantics rather than R's 32-bit integer storage type. Result lists contain diagnostics and are not promised as a long-term serialization format.

Building the frame

The frame is a data frame where each row is one stratum of your target population. It summarizes what you know about each subgroup before sampling. A typical workflow:

Identify strata from a census or register (e.g. provinces, urban/rural areas, age groups).
Look up N: the population count per stratum.
Estimate sd: the standard deviation of your key variable within each stratum (from a pilot survey, a previous census, or expert judgment). If unknown, set sd = 1 everywhere for proportional allocation.
Add mean or p if you want to solve for a target CV.

A minimal frame:

frame <- data.frame(
  stratum = c("Urban", "Rural"),
  N       = c(50000, 120000),
  sd      = c(12, 20)
)

When a design stratifies by several variables (e.g. region $\times$ urbanicity), cross them into one variable first:

frame$stratum <- interaction(frame$region, frame$urban, drop = TRUE)

This ensures that each row maps to exactly one population cell and that the allocation formulas apply to the correct per-stratum N and sd pairs.

Cluster designs within strata

Adding a icc_psu column to the frame turns the allocation into a stratified two-stage design (e.g. enumeration areas then households within each stratum). Under the cluster variance model the problem reduces to the element allocation above with the stratum SD inflated to sd * sqrt(var_ratio_psu * (1 + icc_psu * (n_per_psu * resp_rate - 1))) and, when stage costs are given, a per-element cost of cost_psu / n_per_psu + cost_ssu. The clustering penalty is paid on the take that responds, n_per_psu * resp_rate, since a unit that does not respond contributes no within-cluster observation. The cost is paid on the gross take, since it is issued either way. The two coincide at resp_rate = 1. The whole-cluster operational search reads the same responding take, so a design it accepts is one the continuous reduction also accepts. All solve modes, allocation methods, and constraints work unchanged. The n, cv, and budget modes keep their meanings.

Cluster-mode columns:

icc_psu (required): within-PSU homogeneity per stratum, e.g. from varcomp() with strata. When PSUs are selected with probability proportional to size, pass prob to varcomp() as well.
var_ratio_psu (optional, default 1): variance ratio per stratum.
n_per_psu (optional): fixes the per-stratum take, in the operational whole-unit design as much as in the continuous one. Any NA entries are replaced by the cost-optimal take sqrt(cost_psu / cost_ssu * (1 - icc_psu) / icc_psu), or its stage-corrected form under fpc = "stage" (see "Which correction the variance carries"). Holding the take is what n mode gives up its exact total for. The stratum is fielded as a whole number of clusters of that size, so its element count lands on the nearest multiple.
cost_psu, cost_ssu (together): per-PSU and per-element costs. These are required for budget mode or when n_per_psu is not fixed. They replace unit_cost, which is not allowed in this mode.
N_psu (optional): the number of PSUs available in the stratum.

What `N_psu` does, and what it does not

N_psu is a feasibility constraint and nothing more. It caps the allocation at N_psu * n_per_psu ultimate units and caps the whole-unit design at N_psu clusters, matching the bound the fixed-take path already applies. Leaving it out preserves the unbounded behavior, in which the allocation may ask for more PSUs than a stratum contains.

A take the frame leaves free is not capped in elements under alloc = "optimal" in the n mode, and in the cv and budget modes without domains, where the search moves the take. Once the allocation passes N_psu PSUs at the cost-optimal take, the stratum fields all of them and its take rises to n / N_psu, with variance and cost read at that take. The allocation rule is unchanged, since on both sides of the switch the variance a unit of cost buys is the same, exactly so without a finite population correction. prec_alloc() assesses a free take the same way, so n_alloc(prec_alloc(fit)) reproduces such a design.

The search keeps a free take within the average PSU size N / N_psu, since no PSU fields more units than it holds, even where the rule's continuous take is larger. Without N_psu the take stays at the rule's whole take.

Under fpc = "none" and "stage" the continuous take is the exact cost-optimal take of the variance the allocation reports. Under the default "unit" it is the take without a correction, since there the optimal take moves with the budget. A slightly different take can then buy a little more precision for the same cost, and the whole design, which the search chooses under the correction, can do better than the continuous one.

On its own it does not activate a first-stage finite population correction. Under the default fpc = "unit" no factor keyed to the PSU sampling fraction enters the variance, however large a share of the PSU universe the design takes, and print() says so once that share is appreciable. Supplying N_psu bounds the design. Asking for the correction as well is what fpc = "stage" does.

So at n_psu == N_psu the default still carries between-PSU variance even though every PSU has been selected. That is deliberately conservative and is then not a literal variance. For the same reason, a target this path reports as unreachable at the PSU bound may be reachable under a finite-population first stage, and the error says so rather than claiming the precision is impossible. Under fpc = "stage" the between-PSU component does reach zero there, and the whole variance reaches zero only when every stage is enumerated.

Because taking every PSU leaves the within-PSU take in force, it does not enumerate a stratum, and take_all is refused in cluster mode for the same reason the fixed-take path refuses it.

The whole-unit design keeps its PSU count within N_psu as well. When the strata at their N_psu limits cannot hold n at their cost-optimal whole takes, a take is raised above its continuous value. When no whole take brings the design within half a take of n, n_alloc() stops with an error giving the largest total it can field.

Which correction the variance carries

fpc chooses, and the three choices are nested rather than arbitrary. Write $f_1$ for the PSU sampling fraction and $f_2$ for the within-PSU one, both measured on what responds. Under the equal-size, equal-take variance-component model this path assumes, the two-stage without-replacement variance corrects each component by its own stage:

$$V = (1-f_1)\frac{S_1^2}{a} + (1-f_2)\frac{S_2^2}{am}.$$

"unit" (default) applies the ultimate-unit 1 - n / N to the whole clustering-inflated variance, so the between-PSU component is attenuated by it too. Since $n/N = f_1f_2$, that factor is at least as large as either exact one, so this overstates both components rather than either being understated. It is conservative against the exact without-replacement variance at every sampling fraction.
"stage" applies the expression above, and its three-stage analogue. It needs N_psu, and N_ssu as well at three stages. It is exact for that model and not beyond it: unequal PSU sizes, a size-proportional first stage, and stochastic rather than expected nonresponse each need a more general variance expression, so it is refused with a psu register.
"none" applies no correction at any stage, which is n_cluster()'s default model exactly. Use it to compare the two interfaces at a sampling fraction where the correction would otherwise separate them.

So none $\ge$ unit $\ge$ stage at every sampling fraction, and the three coincide as the fractions vanish.

The default stays "unit" deliberately. "stage" is exact under that model, but it buys the exactness by depending on the mean PSU size $M = N/N_{psu}$, which the with-replacement form does not use at all. Under unequal PSU sizes that mean is a stand-in, and planners usually know it less well than they know icc_psu and the take. Choose "stage" when the PSU sizes really are close to equal and the fractions are large enough to matter.

Under "stage" a free take is the cost-optimal take of the corrected variance. With $r$ the unit response rate it is

$$m = \sqrt{\frac{c_1(1-\delta)}{c_2\,r\,\{\delta - (1-\delta)/M\}}}.$$

The within-PSU fraction lowers the between-PSU term to $\delta - (1-\delta)/M$. When $\delta \le 1/(M+1)$ that term is not positive, a further unit inside a PSU always buys more precision than a further PSU, and the take is the whole PSU, $M$. The take moves well away from the uncorrected one only for small PSUs. At $M = 50$ and $\delta = 0.02$, keeping the uncorrected take would cost about a fifth of the variance at the same budget, and at a few hundred units per PSU the two takes give the same variance to within half a percent.

Note that the two fractions are controlled by different quantities. $f_1$ falls as the design takes a smaller share of the PSU universe, but $f_2$ is the realized take over $M$ and does not move with N at all. A large population reached through a few large PSUs leaves $f_2$ wherever it was, so "stage" is not in general a small correction to "none".

An element allocation has no stages to choose between, so it takes "unit", the correction 1 - n / N, or "none", none at all. Under "none" a stratum variance derived from p is p (1 - p), as in n_prop() with N = Inf, rather than the finite-population N p (1 - p) / (N - 1). The stratum weights still read N (times eligibility_rate), and the caps, take-all strata and base weights stay in frame units, so only the correction is dropped. A take-all stratum then keeps a sampling variance, as a census of PSUs does in n_cluster(). Use "none" to reproduce a plan built on published design effects, which carry the clustering and leave no element correction to apply.

Because icc_psu already accounts for the clustering, leave deff at 1 unless it captures a different source of design effect (e.g. weighting loss). A clustering deff on top of icc_psu would double-count. The constraints min_n_stratum and max_weight stay in element units. For fielding, use the whole-unit design in n_psu_int and n_per_psu_int (n_int = n_psu_int * n_per_psu_int). Its actual field cost and precision are reported in $operational and, in budget mode, never exceed the budget.

Certainty PSUs from a register

N_psu gives only the number of PSUs in a stratum. Supply psu to give their sizes instead. It has one row per PSU with stratum and N, and the sizes must sum to frame$N in every stratum. This two-stage joint-allocation mode requires n_per_psu in frame and icc_psu in measures.

A remainder PSU is drawn for the take, so a PSU holding fewer units than n_per_psu is refused, named by its psu_id or its row. Merge it with a neighboring PSU before planning, as is usual when PSUs are formed (Valliant, Dever and Kreuter 2018, sec. 10.3) and as merge_psus() does, or flag it in certainty, where it is subsampled at the stratum's sampling fraction. A psu_id must name each PSU once.

icc_psu describes the remainder, whose PSUs are selected with probability proportional to size. It is the size-weighted share of the stratum variance that lies between PSUs, which varcomp() returns when prob is each PSU's share of its stratum's units, as in the last example. A stratum with a single PSU has no variance between PSUs, so its icc_psu is taken as 0 whatever is supplied and may be NA. Every other stratum needs a value in $[0, 1]$.

With take $b_h$, sampling fraction $f_h = n_h / N_h$ and certainty_cutoff $c_h$, a PSU is certain when its inclusion probability reaches $c_h$, that is when its size reaches

$$N_{hi} \ge c_h b_h / f_h.$$

The remainder is fielded as a whole number of PSUs drawn with probability proportional to size, and rounding that number up raises every remainder PSU's inclusion probability. A PSU whose probability reaches $c_h$ under that draw is certainty too, and the draw is recomputed without it until no PSU reaches $c_h$, the iterative rule for selection proportional to size. A PSU below the threshold can be certainty this way, and so can a stratum's only PSU. A remainder always draws at least one PSU, even when the rounded certainty takes use up the stratum's allocation, so no part of a stratum is left out of the sample. The operational design is therefore one any fixed-size proportional-to-size method can draw as it stands.

The default $c_h = 1$ takes with certainty only the PSUs the design cannot sample less often. A cutoff below one, such as the 0.80 Valliant, Dever and Kreuter (2018, sec. 9.6) mention, takes more PSUs with certainty and removes their clustering from the design effect, so the sample usually falls. The classification is discrete, so the plan moves in steps as the cutoff changes and need not fall with it. Price the choice with predict.svyplan over certainty_cutoff.

n_psu_per_zone = 2 cuts each stratum's remainder, sorted by psu$zone_order or by size, into zones of about equal total size and draws two PSUs from each, the deep stratification of area samples (Valliant, Dever and Kreuter 2018, sec. 9.5 and Example 3.13). The remainder draw rounds up to a whole number of zones, which can add one PSU in a stratum. The rule above then applies inside each zone. A PSU is certain when its probability within its zone reaches $c_h$, and so is every PSU of a zone, or of a remainder, that holds no more PSUs than it draws. A stratum that draws from its remainder then draws at least two PSUs, so its variance can be estimated. The design effect does not credit the zones, so the planned precision stays conservative.

n_psu_per_zone = 1 draws one PSU from each zone, the finest stratification of the remainder. A zone holding a single PSU is a census, so more PSUs are certainty and the sample tends to be the smallest of the three settings. No zone yields two PSUs, so the variance needs zones collapsed, and $psu$.pair fixes the groups before selection. Zones are paired in zone order within a stratum, with three together where the count is odd. Strata whose remainder is a single zone are grouped with each other in frame row order, so order the rows to put similar strata next to each other. One draw per zone has the with-replacement variance, so the groups are analyzed without a finite population correction. In that form the collapsed estimator overestimates the variance, by more as the zone totals in a group differ (Valliant, Dever and Kreuter 2018, sec. 15.5.3). A without-replacement correction such as Brewer's can turn it into an underestimate when the zones in a group are alike.

The certainty part has no first-stage sampling variance. The remainder has the usual clustering component, and the two are combined into the anticipated design effect. psu$certainty can add a PSU to the certainty part. It cannot remove a PSU either rule makes certainty.

The continuous allocation uses an aggregate planning approximation. The noncertainty part uses a with-replacement clustering approximation, a standard approach in PPS sample-size planning. Each part receives its proportional share of the stratum allocation, so with $N_h^C$ units in certainty PSUs, $N_h^R$ in the remainder, $\delta_h$ for icc_psu and response rate $r_h$, the design effect is

$$d_h = \frac{N_h^C + N_h^R \{1 + \delta_h (b_h r_h - 1)\}}{N_h}.$$

It depends on which PSUs are certainty, not on the allocation. A common stratum-level finite population correction is then applied. The clustering reads the responding take $b_h r_h$, as the N_psu form does, while the threshold and the draw read the issued take $b_h$. A PSU-level response rate and a var_ratio_psu other than 1 have no place in this design effect and are refused.

The field design is chosen whole. For each stratum the planner lists every design these rules give as the sampling fraction runs, and under stage costs also those holding the stratum's largest PSUs in certainty below the threshold, which cuts clustering at the price of their visits. Under unit cost a visit has no price, so no PSU is held below the threshold. Give cost_psu and cost_ssu where holding a large PSU in certainty may be worth its visit. A search then picks one design per stratum, the cheapest meeting every target, and proves it so, or stops after control$max_nodes nodes, 5,000 by default, with the best design found and a lower bound on the cheapest, which $optimization$search reports. With objective and budget it picks the design with the smallest weighted relative variance among those meeting every target in targets whose field cost, visits included, stays within the budget. A budget below the cheapest field design, or below the cheapest one meeting the targets, is an error that gives that design's cost. The continuous allocation in $detail$n solves the aggregate model above under the planned classification, as a relaxation of the field design. $detail gives the split and threshold by stratum, while $psu gives each PSU's classification, its source, and n_take, the take the operational design fields in it. A certainty PSU carries its own whole take at the stratum rate, capped at its size, and a remainder PSU carries n_per_psu. The takes are the operational design's own numbers, so sum(n_take[certainty]) = n_certain_int in every stratum.

A fit solved with a register keeps it in $params$psu. Its presence is the stable test for a certainty-aware fit, and code can key on it rather than on $detail column names.

The operational design is fieldable. Certainty PSUs use their whole take at the stratum rate and the remainder uses whole PSUs at n_per_psu. Thus n_int = n_certain_int + n_psu_draw * n_per_psu in every stratum. $operational reports its count, cost and precision. The rounded certainty takes $q_{hi}$ run above the stratum rate and the remainder gets what is left, so the field design is assessed on its own takes and draw $a_h$, each part a stratum of its own:

$$V_h = S_h^2 \Big[\sum_{i \in C} N_{hi}^2 \Big(\frac{1}{r_h q_{hi}} - \frac{1}{N_{hi}}\Big) + (N_h^R)^2 \{1 + \delta_h (b_h r_h - 1)\} \Big(\frac{1}{r_h a_h b_h} - \frac{1}{N_h^R}\Big)\Big].$$

At the stratum rate the expression equals the aggregate one above. It is the variance the planner sizes, so every target holds on the design that is fielded. A census of every PSU gives each target its smallest variance, so when even that misses a target, n_alloc() stops with an error that says so.

Use cost_psu and cost_ssu together to price PSU visits and interviews. A certainty PSU's visit is paid whatever its take, so the solver prices an interview at cost_ssu + (N_h^R / N_h) * cost_psu / n_per_psu, and $params$achieved$cost and $optimization$cost add the certainty visits as a fixed cost, as do the stratum costs summary() reports. The optimality certificate, its gap and its multipliers refer to the marginal-cost problem without that constant, which has the same optimum, and hold for the planned classification, not across classifications. In budget mode the continuous solve has the budget less the certainty visits. The field design prices its visits exactly, so near the cheapest design meeting the targets the aggregate model can need more than the budget. It is then solved at the cost of its own cheapest design meeting them, and $optimization$budget_residual is negative. predict.svyplan compares values of n_per_psu, certainty_cutoff, n_psu_per_zone and, in budget mode, budget, and reports each row's search status.

Domains vs. strata

Domains are specified via the domains parameter. Domain columns partition strata into sub-populations. Each domain groups one or more strata. When cv is specified, every domain must meet the target, and since a domain's precision depends on its own strata only, the smallest design that does so is each domain at its own smallest total. Each domain is allocated by the rule inside itself, as if it were the frame, and its total is the least that meets the target there. The result is the minimum joint allocation reaches with a cv target in every domain. National precision is reported in $cv and is not a target. To require it as well, state national and domain targets together in joint allocation (measures and targets). Under alloc = "power", alloc_q shapes the shares inside a domain, and the domain totals come from their targets.

In n or budget mode, domains affect reporting only. Per-domain precision metrics appear in $domains, but the allocation itself treats all strata globally.

min_n_domain adds a floor on each domain's total, in every mode. A domain whose share of the allocation falls short of its floor is held at the floor and allocated inside itself by the same rule, and the rest of the total is spread over the other domains. This repeats until no domain falls short. Holding one domain lowers the others' shares, so a domain already above its floor can end up held too. The result is the best allocation by the rule under the floors. In cv and budget mode the floors raise the total or the cost only where they bind. The whole-unit design keeps every floor, and in n and budget mode it takes the units it adds to a short domain from domains with room. With n_step it adds and removes whole steps, and in n mode it reaches the total by moving several strata together when no single step fits the excess, as 36 and 36 for steps of 12, a total of 72 and floors of 25. $domains then reports each floor as .min_n and whether it binds as .floor_binding.

Allocation methods

Allocation is controlled by the alloc parameter (same methods as strata_bound()):

proportional: $n_h \propto N_h / r_h$
neyman: $n_h \propto N_h S_h \sqrt{d_h / r_h}$
optimal: $n_h \propto N_h S_h \sqrt{d_h / r_h} / \sqrt{c_h}$
power: Bankier (1988) allocation, $n_h \propto (S_h/|\mu_h|) X_h^q\sqrt{d_h/r_h}$, where $X_h$ is frame$alloc_measure, defaulting to $N_h|\mu_h|$. Every stratum needs a finite nonzero mean or p. The population CV $S_h/|\mu_h|$ is distinct from the target cv, which describes the relative standard error of an estimate. Zero-variance strata have zero allocation weight and remain subject to the allocation bounds. If all variances are zero, the existing fallback allocates in proportion to population size, with a warning. Use explicit domain targets when comparable domain precision is required.

The design effect $d_h$ and response rate $r_h$ enter only when they vary by stratum. A value shared by every stratum is a constant factor and cancels out of a proportional weighting, leaving the classical rules above it. Neyman and optimal allocation carry $\sqrt{d_h/r_h}$ from minimizing $\sum W_h^2S_h^2d_h/(r_hn_h)$ under a constraint on issued count or cost, respectively. The power rule retains that adjustment while weighting population CVs by alloc_measure^alloc_q. Proportional carries $1/r_h$ instead, and no $d_h$ at all, because it is a count rule, not a variance optimum, and its purpose is a self-weighting sample, so it is the responding sample it holds proportional to $N_h$.

Stratum allocations are rounded to integers using the ORIC method (Cont and Heidari, 2015). Constraints (min_n_stratum, max_weight, take_all) are enforced via recursive Neyman allocation (RNA, Wesolowski et al., 2021).

When budget is specified, the algorithm finds the maximum affordable allocation under unit costs.

Sample sizes in steps

A frame column n_step makes each stratum's n_int a multiple of its step. The continuous allocation is solved within bounds moved to whole steps, then rounded in steps. In n mode each stratum takes the step nearest its share, and then moves by at most one step, jointly across strata, to reach n exactly. A total the steps cannot reach within that window is an error naming the nearest totals they do reach, so n = 550 with steps of 4 and 8 names 548 and 552. A cv target rounds each stratum up to a step. A budget adds and removes whole steps. The searched modes search whole steps. With one step for every stratum the result is the allocation in sets of that size multiplied back. prec_alloc() assesses any n on such a frame, so a previous round's allocation can be compared with one in steps.

References

Valliant, R., Dever, J. A., & Kreuter, F. (2018). Practical Tools for Designing and Weighting Survey Samples (2nd ed.). Springer. Chapters 5 and 9.

Bankier, M. D. (1988). Power allocations: determining sample sizes for subnational areas. The American Statistician, 42(3), 174–177.

Bethel, J. (1989). Sample allocation in multivariate surveys. Survey Methodology, 15(1), 47–57.

Cont, R. and Heidari, M. (2015). Optimal rounding under integer constraints. arXiv preprint arXiv:1501.00014.

Wesolowski, J., Wieczorkowski, R. and Wojciak, W. (2021). Optimality of the recursive Neyman allocation. Journal of Survey Statistics and Methodology, 10(5), 1263–1275.

Examples

frame <- data.frame(
  stratum = c("A", "B", "C"),
  N    = c(4000, 3000, 3000),
  sd   = c(10, 15, 8),
  mean = c(50, 60, 55),
  unit_cost = c(1, 1.5, 1)
)

n_alloc(frame, n = 600)
#> Stratum allocation (neyman, 3 strata)
#> field design: n = 600, cv = 0.0079, cost = 724
#> continuous design: n = 600, cv = 0.0079, se = 0.4305
#> (deff = 1)
#> design df = 597
n_alloc(frame, cv = 0.03)
#> Stratum allocation (neyman, 3 strata)
#> field design: n = 45, cv = 0.0298, cost = 54
#> continuous design: n = 44.23479, cv = 0.0300, se = 1.6350
#> (deff = 1)
#> design df = 42

frame_constraints <- transform(
  frame,
  max_weight = c(25, 20, NA),
  take_all = c(FALSE, FALSE, TRUE)
)

n_alloc(frame_constraints, budget = 3500, alloc = "optimal", min_n_stratum = 40)
#> Stratum allocation (optimal, 3 strata)
#> field design: n = 3403, cv = 0.0076, cost = 3500
#> continuous design: n = 3403.425, cv = 0.0076, se = 0.4125
#> (min_n_stratum = 40, deff = 1)
#> design df = 401

frame_domains <- data.frame(
  province = c("North", "North", "South", "South"),
  stratum = c("Urban", "Rural", "Urban", "Rural"),
  N    = c(2000, 3000, 1800, 3200),
  sd   = c(12, 18, 10, 16),
  mean = c(55, 48, 58, 50)
)

n_alloc(frame_domains, domains = "province",
       cv = 0.04, alloc = "power", alloc_q = 0.3)
#> Stratum allocation (power, 4 strata)
#> field design: n = 104, cv = 0.0280, cost = 104
#> continuous design: n = 101.8385, cv = 0.0283, se = 1.4666
#> (deff = 1)
#> design df = 100
#> Domains: 2
#> 
#>  province .domain .n       .se    .moe     .rmoe      .cv    .cost
#>  North    5_North 58.67031 2.0320 3.982647 0.07839856 0.0400 59   
#>  South    5_South 43.16817 2.1152 4.145716 0.07839856 0.0400 43   

# A publication floor of 480 per province, which holds the South, and
# one of 550 for the South only
n_alloc(frame_domains, domains = "province", n = 1000,
        min_n_domain = 480)$domains
#>   province .domain  .n       .se     .moe      .rmoe        .cv .cost .min_n
#> 1    North 5_North 520 0.6462198 1.266567 0.02493243 0.01272086   520    480
#> 2    South 5_South 480 0.5992373 1.174484 0.02221035 0.01133202   480    480
#>   .floor_binding
#> 1          FALSE
#> 2           TRUE
south <- data.frame(domain = "province", level = "South", min_n = 550)
n_alloc(frame_domains, domains = "province", n = 1000,
        min_n_domain = south)$domains
#>   province .domain  .n       .se     .moe     .rmoe        .cv .cost .min_n
#> 1    North 5_North 450 0.7002857 1.372535 0.0270184 0.01378515   450     NA
#> 2    South 5_South 550 0.5552447 1.088260 0.0205798 0.01050009   550    550
#>   .floor_binding
#> 1          FALSE
#> 2           TRUE

# Power allocation on a size measure alone, as DHS and MICS apply it
sizes <- data.frame(stratum = c("A", "B", "C"), N = c(4e5, 1e5, 2.5e4))
n_alloc(sizes, n = 60, alloc = "power", alloc_q = 0.25, min_n_stratum = 15)
#> Stratum allocation (power, 3 strata)
#> field design: n = 60, cost = 60
#> continuous design: n = 60
#> (min_n_stratum = 15, deff = 1)
#> precision: not assessed (no 'sd', 'var' or 'p' in frame)
#> design df = 57

# A floor for one stratum only, such as the clusters a small domain needs
# for its own estimates. The others share the rest by the power rule.
n_alloc(sizes, n = 60, alloc = "power", alloc_q = 0.25,
        min_n_stratum = c(NA, NA, 20))
#> Stratum allocation (power, 3 strata)
#> field design: n = 60, cost = 60
#> continuous design: n = 60
#> (min_n_stratum per stratum, deff = 1)
#> precision: not assessed (no 'sd', 'var' or 'p' in frame)
#> design df = 57

# The minimum is a bound. To give each stratum its minimum and share the
# rest by size on top of it, allocate what is left and add the minimums.
floors <- c(10, 10, 10)
n_alloc(sizes, n = 60 - sum(floors), alloc = "proportional")$detail$n +
  floors
#> [1] 32.85714 15.71429 11.42857

# PSU counts in steps: 4 rotation groups, and 8 in strata split in two
psus <- data.frame(
  stratum = c("North", "South", "East", "West"),
  N = c(120000, 95000, 60000, 30000),
  min_n_stratum = c(16, 12, 12, 13),
  n_step = c(4, 8, 4, 4)
)
n_alloc(psus, n = 120, alloc = "proportional")$detail[
  c("stratum", "n", "n_int", "n_step")]
#>   stratum        n n_int n_step
#> 1   North 45.38182    48      4
#> 2   South 35.92727    32      8
#> 3    East 22.69091    24      4
#> 4    West 16.00000    16      4

# A proportion supplies its own variance
rates <- data.frame(stratum = c("A", "B", "C"), N = c(4e5, 1e5, 2.5e4),
                    p = c(0.48, 0.60, 0.16))
n_alloc(rates, cv = 0.05)
#> Stratum allocation (neyman, 3 strata)
#> field design: n = 406, cv = 0.0500, cost = 406
#> continuous design: n = 405.9453, cv = 0.0500, se = 0.0244
#> (deff = 1)
#> design df = 403

# Joint targets for two indicators and an overlapping domain
joint_frame <- data.frame(
  stratum = c("North urban", "North rural", "South urban"),
  region = c("North", "North", "South"),
  N = c(1000, 1800, 1200)
)
joint_measures <- data.frame(
  stratum = rep(joint_frame$stratum, 2),
  name = rep(c("coverage", "income"), each = 3),
  p = c(0.5, 0.4, 0.6, NA, NA, NA),
  mean = c(NA, NA, NA, 50, 55, 60),
  sd = c(NA, NA, NA, 10, 12, 15)
)
joint_targets <- data.frame(
  name = c("coverage", "income"),
  domain = c(".overall", "region"),
  level = c(NA, "North"),
  cv = c(0.06, NA),
  moe = c(NA, 3)
)
joint <- n_alloc(
  joint_frame, measures = joint_measures, targets = joint_targets
)
prec_alloc(joint, n = joint$detail$n_int)
#> Joint allocation precision (2 constraints)
#> targets: all pass

# The best design a fixed budget can buy, subject to the same targets
bought <- n_alloc(
  joint_frame, measures = joint_measures, targets = joint_targets,
  objective = "income", budget = 4000
)
bought$objective
#>         component   name   domain level priority .relvar .cv .contribution
#> 1 income@.overall income .overall  <NA>        1       0   0             0
#>   .share
#> 1     NA

# Unequal priorities, and an objective on a specific domain
priorities <- data.frame(
  name     = c("income", "coverage"),
  domain   = c(".overall", "region"),
  level    = c(NA, "North"),
  priority = c(2, 1)
)
n_alloc(
  joint_frame, measures = joint_measures, targets = joint_targets,
  objective = priorities, budget = 4000
)
#> Joint constrained allocation (Bethel)
#> question: best design affordable within a budget of 4000
#> field design: n = 4000, cost = 4000 (2 targets, all pass)
#> continuous optimum: n = 4000, cost = 4000
#> objective: weighted relative variance 7.63819e-20 continuous, 7.63819e-20 operational
#> the budget is not binding: the allocation sits at its upper bounds
#> 
#>  component             priority .cv          .share
#>  income@.overall       2        0.000000e+00 0     
#>  coverage@region=North 1        2.763728e-10 1     
#> 
#> active bounds: 0 lower, 3 upper

# A domain without a level column targets every level of it
every_region <- data.frame(
  name = c("coverage", "income"),
  domain = c(".overall", "region"),
  cv = c(0.06, 0.05)
)
by_region <- n_alloc(
  joint_frame, measures = joint_measures, targets = every_region
)
by_region$constraints[, c("constraint", ".target", ".achieved")]
#>               constraint .target  .achieved
#> 1   coverage@.overall:cv    0.06 0.06000000
#> 2 income@region=North:cv    0.05 0.01502439
#> 3 income@region=South:cv    0.05 0.02703700

# What would other budgets buy? (cost-versus-objective frontier)
predict(bought, data.frame(budget = c(3000, 4000, 5000)))
#>   budget    n cost objective_value          cv n_int cost_int .binding
#> 1   3000 3000 3000    3.907646e-06 0.001976777  3000     3000     TRUE
#> 2   4000 4000 4000    0.000000e+00 0.000000000  4000     4000    FALSE
#> 3   5000 4000 4000    0.000000e+00 0.000000000  4000     4000    FALSE
#>   .feasible
#> 1      TRUE
#> 2      TRUE
#> 3      TRUE

# The same target tables with fixed two-stage takes
joint_frame$N_psu <- c(100, 150, 100)
joint_frame$n_per_psu <- c(8, 10, 12)
joint_frame$cost_psu <- c(300, 400, 450)
joint_frame$cost_ssu <- c(25, 30, 35)
joint_measures$icc_psu <- rep(c(0.03, 0.05, 0.08), 2)
joint_measures$var_ratio_psu <- 1
joint_cluster <- n_alloc(
  joint_frame, measures = joint_measures, targets = joint_targets
)
joint_cluster$detail[, c("stratum", "n_psu_int", "n_per_psu", "n_int")]
#>       stratum n_psu_int n_per_psu n_int
#> 1 North urban        11         8    88
#> 2 North rural        17        10   170
#> 3 South urban        11        12   132

# Each indicator measured on its own eligible persons per household
joint_measures$eligible_per_unit <- rep(c(0.41, 2.6), each = 3)
n_alloc(joint_frame, measures = joint_measures, targets = joint_targets)
#> Joint constrained allocation (Bethel)
#> question: cheapest design meeting every precision target
#> field design: n = 684, cost = 47160 (2 targets, all pass)
#> continuous optimum: n = 683.0394, cost = 47038 (integerizing costs +0.26%)
#> binding: coverage@.overall:cv (target 0.06, achieved 0.06)

# Fixed three-stage takes: PSUs, SSUs within PSU, ultimate units within SSU
joint_frame$N_ssu <- c(600, 1000, 750) # aggregate SSUs by stratum
joint_frame$n_per_ssu <- c(3, 4, 3)
joint_frame$cost_tsu <- c(5, 6, 7)
joint_measures$icc_ssu <- rep(c(0.10, 0.08, 0.12), 2)
joint_measures$var_ratio_ssu <- 1
joint_three <- n_alloc(
  joint_frame, measures = joint_measures, targets = joint_targets,
  min_n_stratum = 20
)
joint_three$detail[, c(
  "stratum", "n_psu_int", "n_per_psu", "n_per_ssu", "n_int"
)]
#>       stratum n_psu_int n_per_psu n_per_ssu n_int
#> 1 North urban         9         8         3   216
#> 2 North rural        12        10         4   480
#> 3 South urban         8        12         3   288
prec_alloc(joint_three, n = joint_three$detail$n_int)
#> Joint allocation precision (2 constraints)
#> targets: all pass

# Stratified two-stage design (EAs then households per stratum)
frame_cluster <- data.frame(
  stratum   = c("Urban", "Rural"),
  N         = c(50000, 150000),
  sd        = c(0.45, 0.48),
  mean      = c(0.35, 0.25),
  icc_psu = c(0.03, 0.08),
  cost_psu  = c(300, 600),
  cost_ssu  = c(40, 60)
)

n_alloc(frame_cluster, cv = 0.05)
#> Stratum allocation (neyman, two-stage, 2 strata)
#> field design: n = 2005, n_psu = 170, cv = 0.0500, cost = 205560
#> continuous design: n = 1979.882, cv = 0.0500, se = 0.0138
#> (deff = 1)
#> design df = 168

# Small PSUs of 49 units each: the stage correction moves the
# cost-optimal take from about 6.3 to 7
plots <- data.frame(N = 4900, N_psu = 100, sd = sqrt(54), mean = 10,
                    icc_psu = 5 / 54, cost_psu = 12, cost_ssu = 3)
n_alloc(plots, budget = 1000, fpc = "none")$detail$n_per_psu
#> [1] 6.26099
n_alloc(plots, budget = 1000, fpc = "stage")$detail$n_per_psu
#> [1] 7

# Certainty PSUs from a register that lists every PSU with its size.
# The "City" stratum is a single PSU.
set.seed(1)
size <- c(round(rlnorm(60, log(150), 0.8)), 4000)
register <- data.frame(
  psu_id = seq_along(size),
  stratum = rep(c("North", "South", "City"), c(30, 30, 1)),
  N = size
)
pop <- data.frame(
  stratum = rep(register$stratum, size),
  psu_id = rep(register$psu_id, size),
  y = rbinom(sum(size), 1, rep(plogis(rnorm(61, -0.5, 0.4)), size))
)
# icc_psu under size-proportional selection: each PSU's share of its stratum
pop$share <- ave(pop$y, pop$psu_id, FUN = length) /
  ave(pop$y, pop$stratum, FUN = length)
vc <- varcomp(y ~ psu_id, strata = ~stratum, prob = ~share,
              data = subset(pop, stratum != "City"))

reg_frame <- aggregate(N ~ stratum, register, sum)
reg_frame$n_per_psu <- 10
reg_measures <- aggregate(cbind(p = y) ~ stratum, pop, mean)
reg_measures$name <- "y"
# NA for "City": a single PSU has no between-PSU variance
reg_measures$icc_psu <- vc$strata$icc_psu[
  match(reg_measures$stratum, vc$strata$stratum)
]
reg_measures
#>   stratum         p name    icc_psu
#> 1    City 0.3330000    y         NA
#> 2   North 0.3906302    y 0.03284821
#> 3   South 0.3829961    y 0.03343549

reg_fit <- n_alloc(
  reg_frame, measures = reg_measures,
  targets = data.frame(name = "y", cv = 0.05), psu = register
)
reg_fit$detail[, c("stratum", "n_psu_certain", "n_psu_draw", "n_int")]
#>   stratum n_psu_certain n_psu_draw n_int
#> 1    City             1          0   155
#> 2   North            12         10   272
#> 3   South            10         11   270

# Relaxing the certainty rule to an inclusion probability below one
predict(reg_fit, data.frame(certainty_cutoff = c(1, 0.9, 0.8)))
#>   certainty_cutoff n_psu_certain n_psu_draw        n     cost n_int cost_int
#> 1              1.0            23         21 696.2411 696.2411   697      697
#> 2              0.9            27         17 687.0380 687.0380   690      690
#> 3              0.8            30         14 680.1932 680.1932   682      682
#>   .status .feasible
#> 1 optimal      TRUE
#> 2 optimal      TRUE
#> 3 optimal      TRUE

# No zones, one PSU per zone, and two PSUs per zone
predict(reg_fit, data.frame(n_psu_per_zone = c(NA, 1, 2)))
#>   n_psu_per_zone n_psu_certain n_psu_draw n_zone        n     cost n_int
#> 1             NA            23         21      0 696.2411 696.2411   697
#> 2              1            40          7      7 663.5245 663.5245   667
#> 3              2            34         12      6 673.0974 673.0974   674
#>   cost_int .status .feasible
#> 1      697 optimal      TRUE
#> 2      667 optimal      TRUE
#> 3      674 optimal      TRUE

# The register design a field budget buys, with visits and interviews
# priced and the cv of 0.05 kept as a hard target
reg_frame$cost_psu <- 300
reg_frame$cost_ssu <- 20
reg_buy <- n_alloc(
  reg_frame, measures = reg_measures,
  targets = data.frame(name = "y", cv = 0.05), objective = "y",
  budget = 30000, psu = register
)
reg_buy$operational[c("cost", "objective_value", "budget_residual")]
#> $cost
#> [1] 30000
#> 
#> $objective_value
#> [1] 0.002065066
#> 
#> $budget_residual
#> [1] 0
#>