Package 'rmoriebricklayer'

Description

Works out how many admissions a given average daily population and average length of stay imply over a period. Published statistics often give the population held and the average stay but not the number admitted (or the other way round); this recovers the missing quantity from the identity that links all three.

Usage

admits(adp, alos, t = 365)

Arguments

Details

Rearranging the identity $\bar{X}_{hc} = N_a \bar{X}_t / t$ (Lakner 1976, p.18-20) for $N_a$. Useful when two of the three quantities are published and the third is not. The identity holds exactly only in a steady state (admissions and stays roughly constant across the period); with a population that is growing or shrinking quickly the result is an approximation.

Errors when any argument is not numeric, is NA or Inf, is not positive, or has more than one value.

Value

A numeric scalar: the implied number of admissions over the period, $\mathrm{adp} \times t / \mathrm{alos}$. It is not rounded. It has class c("rmbl_admits", "rmbl_measure", "numeric") and an attribute "lakner" (adp, alos and t), so it prints as a box with its inputs, formula and caveats; arithmetic, round(), sum(), mean(), [, c() and as.numeric() give a plain number.

See Also

adp(), alos(), stkflow()

Examples

# Lakner (p.20) runs this the other way: t = 365, an average daily
# population of 25 and 1,750 admissions imply a stay of 5.2 days.
alosimplied <- 25 * 365 / 1750
round(alosimplied, 1)

# and back again
admits(25, alosimplied)

# a 90-day quarter with 40 held on average and a 12-day average stay
admits(40, 12, t = 90)

# an average stay of zero is refused
try(admits(25, 0))

Description

The mean number of person-days served per day over a period: a STOCK, answering how many people are held at one time. Use it to report how much confinement (or bed occupancy, shelter use, open caseload) a system used over a year or any other window, from the days each person spent in it. It is the measure capacity planning needs, because it counts beds in use rather than people who passed through.

Usage

adp(days, t = 365)

Arguments

Details

This is Lakner's $\bar{X}_{hc} = \sum X_i / t$ (1976, p.15), where $X_i$ is the days served by person $i$ inside the period. The denominator is TIME, which is what makes it a stock. Compare alos(), whose denominator is people: the two are linked by $\mathrm{ADP} = N \times \mathrm{ALOS} / t$ (Little's law), so admits() and stkflow() can move between them.

Only days served inside the period should be counted: a stay that began before the period or ends after it contributes its days within the window.

Errors when days is not numeric (text and factors are refused), is empty, has NA, Inf or a negative value, or when t is not a single positive number.

Value

A numeric scalar: person-days per day, that is the average number of people held on any one day of the period. It has class c("rmbl_adp", "rmbl_measure", "numeric") and an attribute "lakner" (the total days and t), so it prints as a box with its inputs, formula and caveats; arithmetic, round(), sum(), mean(), [, c() and as.numeric() give a plain number.

References

Lakner, E. (1976) A Manual of Statistical Sampling Methods for Corrections Planners. University of Illinois at Urbana-Champaign.

See Also

alos(), stkflow(), adpcnts()

Examples

# Lakner's own worked example (p.15): 3,000 inmates served 13,500
# detention days in a year.
adp(13500)

# per-person days give the same total
adp(c(10, 20, 30), t = 30)

# a leap year, with the period length taken from the dates
adp(13500, t = perdays("2024-01-01", "2024-12-31"))

# negative days are refused
try(adp(c(10, -2)))

Description

When only a count of people present on certain days is available – not a record per person – the person-days over the period are the mean of those counts scaled to the period's length. Use it with adp() or alos() when the source is a daily or weekly count sheet rather than admission and release dates.

Usage

adpcnts(cnts, t = 365)

Arguments

Details

Lakner's $X_t = (\frac{1}{C}\sum N_i) t$ (1976, p.21). The counts need not cover every day: counting on weekdays only is the usual case, and the mean of the days counted stands in for the days not counted. That substitution is an assumption, and it fails if the days counted differ systematically from the days missed – weekday-only counting where weekend admissions are released before Monday, for instance.

Value

A numeric scalar: the estimated person-days over the period, mean(counts) * t. Despite the function's name this is a total of person-days, not an average daily population; divide by t (or pass it to adp()) for the population held on an average day. It has class c("rmbl_adpcnts", "rmbl_measure", "numeric") and an attribute "lakner" (the number of counts, their mean and t), so it prints as a box with its inputs, formula and caveats; arithmetic, round(), sum(), mean(), [, c() and as.numeric() give a plain number, and adp() takes it as days.

See Also

adp()

Examples

# Lakner (p.21): counts on 255 days summing to 34,935 imply just over
# fifty thousand detention days across the year.
adpcnts(rep(34935 / 255, 255))

# weekly counts over a 31-day month, then the average daily population
pd <- adpcnts(c(120, 130, 125, 128), t = 31)
pd
adp(pd, t = 31)

# a missing count must be dropped or filled first
try(adpcnts(c(120, NA, 125)))

Description

Asks a language model to plan a brick-proof, reproducible data bundle, through llmask: your own key, your own endpoint or a local Ollama server, with the packages' help at hand (the model reads the help pages of the functions it proposes before it answers). Use it to get a first draft of the steps and R code for packaging an analysis: finding the portal's download URL, fetching with a Wayback fallback, recording hashes and a manifest, writing a synthetic stand-in, then packaging, signing and verifying the bundle.

Usage

agentbndl(req, model = NULL, route = NULL)

Arguments

Details

The request is sent with a fixed system prompt that describes that task and asks for a short plan and R code. The route is chosen as in llmask(); when no route is set up, the function returns a message saying what to set up (a provider key, a local Ollama or your own endpoint) instead of an error, and makes no network request. When a route is set up, the request and the prompt are sent to that provider or server. The answer is model output: check the code it proposes against the help pages before running it.

Errors: req not a single non-empty string; route not one of the names listed; route = "hosted" (the hosted tier is retired); and errors raised by llmask().

Value

Character scalar: the answer, or a message saying what to set up when no route answers.

Examples

# Needs a language model: a provider key (e.g. ANTHROPIC_API_KEY), your own
# endpoint or a local Ollama.
agentbndl("scaffold a bundle for analysis.R from the Toronto CKAN data")
agentbndl("sign and verify my capsule", route = "ollama")


# With every route switched off (no key in the environment and no saved
# settings: the settings folder is pointed at an empty one) the call returns
# a setup hint, not an error, and reaches no network -- safe to run anywhere:
keys <- c(
  "OLLAMA_HOST",
  "MORIE_LLM_BASE_URL",
  "XDG_CONFIG_HOME",
  llmprvd()$key_env,
  llmprvd()$alt_env[nzchar(llmprvd()$alt_env)]
)
old <- Sys.getenv(keys, unset = NA)
Sys.setenv(
  OLLAMA_HOST = "off",
  MORIE_LLM_BASE_URL = "off",
  XDG_CONFIG_HOME = tempfile()
)
Sys.unsetenv(setdiff(
  keys,
  c("OLLAMA_HOST", "MORIE_LLM_BASE_URL", "XDG_CONFIG_HOME")
))
cat(agentbndl("hello"))
for (v in keys) {
  if (is.na(old[[v]])) {
    Sys.unsetenv(v)
  } else {
    do.call(Sys.setenv, as.list(old[v]))
  }
}

# argument errors come first
try(agentbndl(""))
try(agentbndl("plan a bundle", route = "hosted"))

Description

The rate a population would have with the age structure of a standard population, for one or more groups, with Fay and Feuer's gamma confidence interval.

Usage

agestd(events, pop, std, grp = NULL, per = 1e+05, conflvl = 0.95)

Arguments

Details

With $d_i$ events and $n_i$ people in age group $i$ and standard weights $w_i$ (the standard population's shares), the rate is $y = \sum w_i d_i / n_i$, its variance $v = \sum (w_i / n_i)^2 d_i$, and with $w_M = \max_i w_i / n_i$ the limits are $\frac{v}{2y} \chi^2_{2y^2/v}(\alpha/2)$ and $\frac{v + w_M^2}{2(y + w_M)} \chi^2_{2(y + w_M)^2/(v + w_M^2)}(1 - \alpha/2)$ (Fay and Feuer 1997). The interval is conservative (coverage at least nominal in Fay and Feuer's simulations) and reduces to the exact Garwood interval for a single age group: suited to the few events of small communities. With no events the lower limit is 0. Rates are per per people. Age groups are matched to std by position: order the rows of each group by age, in the order of std.

The choice of standard matters: Statistics Canada's comparisons of First Nations, Metis and Inuit mortality standardise to the Indigenous population's own age structure (Kumar and Tjepkema 2019), which is younger than Canada's.

Value

A data frame with group, events, rate, lower, upper, and the attributes var (the variances, on the per-person scale) for rategap.

References

Fay MP, Feuer EJ (1997). Confidence intervals for directly standardized rates: a method based on the gamma distribution. Statistics in Medicine 16(7):791-801.
DOI 10.1002/(sici)1097-0258(19970415)16:7<791::aid-sim500>3.0.co;2-#

Kumar MB, Tjepkema M (2019). Suicide among First Nations people, Metis and Inuit (2011-2016): Findings from the 2011 Canadian Census Health and Environment Cohort (CanCHEC). Statistics Canada 99-011-X2019001.

See Also

rategap.

Examples

# made-up counts for one small community, by five age groups, to show the
# computation
std <- c(30000, 25000, 20000, 15000, 10000) # a younger standard population
agestd(
  events = c(4, 9, 11, 6, 3),
  pop = c(9000, 8000, 6000, 4000, 2000),
  std = std
)

Description

The average effect of a binary treatment (on everyone, or on the treated) from one row per person, with the outcome and the propensity fitted by learners on other folds, combined in the doubly robust (Neyman-orthogonal) score of Robins, Rotnitzky and Zhao; with survey or sampling weights if given.

Usage

aipwcf(
  data,
  y,
  trt,
  covs,
  estd = c("ATE", "ATT"),
  lrn = "lasso",
  basis = c("linear", "poly2"),
  w = NULL,
  des = NULL,
  folds = 5L,
  trim = 0.01,
  seed = 1L,
  conflvl = 0.95
)

Arguments

Details

With $\hat\mu_d(X)$ the outcome fitted among those with treatment $d$ and $\hat e(X)$ the propensity, both out of fold, the score for the average effect is

$\psi = \hat\mu_1 - \hat\mu_0 + \frac{D (Y - \hat\mu_1)}{\hat e} - \frac{(1 - D)(Y - \hat\mu_0)}{1 - \hat e} - \theta,$

and for the effect on the treated

$\psi = \frac{D (Y - \hat\mu_0) - \hat e (1 - D)(Y - \hat\mu_0)/(1 - \hat e)}{P(D = 1)} - \frac{D}{P(D = 1)}\theta.$

$\hat\theta$ solves the weighted mean of the score at zero; its variance is the weighted mean of $\psi^2$ over $n$, the influence-function variance of Chernozhukov et al. (2018). Propensities are clipped to [trim, 1 - trim]. The estimator is consistent when either the outcome or the propensity model is right, and with cross-fitting the learners can be flexible (lasso or ridge on a second-order basis, or any learner function).

The average effect equals DoubleML's interactive regression model given the same folds and learners; for the effect on the treated, DoubleML estimates $P(D = 1)$ within each fold and this function in the whole sample, so the two differ slightly.

The propensity learner is logistic regression on the covariates (with basis = "poly2", on their squares and products too). With w the point estimate is design-weighted, but the standard error treats the weights as fixed and the sample as independent.

For a complex survey, pass the design from svydes instead: its weight gives the point estimate and its replicate weights (or strata and primary sampling units) the standard error, with the variance formula of the replicate type, the confidence interval on the design's degrees of freedom. The learners are fitted once, without the weights, which enter only the average of the scores, so each replicate estimate is the replicate-weighted average of the same cross-fitted scores: the same numbers as re-running the function on every replicate weight, from one fit instead of one per replicate. Rows with a missing value are left out as a domain (their scores count as zero in every replicate), not dropped from the design.

Value

A data frame of class rmbl_aipw: estimand, estimate, se, lower, upper, p_value; attributes n, treated, psrange, trimmed (how many propensities were clipped), variance ("influence", or the design's type) and propensity (a data frame of each row's out-of-fold propensity before clipping, its treatment and its weight, which plot() draws).

References

Robins JM, Rotnitzky A, Zhao LP (1994). Estimation of regression coefficients when some regressors are not always observed. Journal of the American Statistical Association 89(427):846-866. DOI 10.1080/01621459.1994.10476818

Chernozhukov V et al. (2018), see plndml.

See Also

plrcf, plndml, baltab.

Examples

set.seed(2)
n <- 800
d <- data.frame(age = runif(n, 18, 70), income = rlnorm(n, 10, 0.5))
d$belong <- rbinom(n, 1, plogis(-2 + 0.03 * d$age + 0.2 * log(d$income)))
d$wellbeing <- 50 + 5 * d$belong + 0.2 * d$age + 2 * log(d$income) +
  rnorm(n, 0, 8)
aipwcf(d, "wellbeing", "belong", c("age", "income"), lrn = "ols")
# with 50 bootstrap replicate weights: one fit, the replicate standard error
d$w <- runif(n, 0.5, 2)
reps <- sapply(1:50, function(b) d$w * rpois(n, 1))
colnames(reps) <- paste0("bsw", 1:50)
d <- cbind(d, reps)
des <- svydes(d, w = "w", reps = colnames(reps))
aipwcf(d, "wellbeing", "belong", c("age", "income"), lrn = "ols",
              des = des)

Description

The mean number of person-days served per person: the FLOW side of the same person-days. Use it to say how long a typical admission lasts, and together with adp() to tell whether a change in the population held came from more people or from longer stays. For the full distribution of stays (median, spread, interval) use staysmry().

Usage

alos(days, n)

Arguments

Details

Lakner's $\bar{X}_t = \sum X_i / N'$ (1976, p.16), with a caveat worth repeating (p.16-17): the period must be longer than the longest stay people actually serve, or the average is biased DOWNWARD, since the longest stays are the ones that fail to finish inside the window. For short-stay facilities a year is comfortable; for long sentences it is not, and the period has to be set from the records.

Errors when days is not numeric, is empty, or has NA, Inf or negative values, and when n is not a single positive number.

Value

A numeric scalar: days per person, $\sum X_i / n$. It has class c("rmbl_alos", "rmbl_measure", "numeric") and an attribute "lakner" (the total days and n), so it prints as a box with its inputs, formula and caveats; arithmetic, round(), sum(), mean(), [, c() and as.numeric() give a plain number.

References

Lakner, E. (1976) A Manual of Statistical Sampling Methods for Corrections Planners. University of Illinois at Urbana-Champaign.

See Also

adp(), stkflow()

Examples

# Lakner (p.17): 2,700 inmates admitted and released served 12,150
# detention days between them.
alos(12150, 2700)

# from one value per released person
stays <- c(2, 5, 1, 30, 4, 7)
alos(stays, length(stays))

# zero people cannot have an average stay
try(alos(12150, 0))

Description

The questions asked of a published table are always the same: what is in it, what changed between periods and how sure can one be, what the rates are if there is an exposure, whether a short series trends, and whether this release differs from the one held before. Each has a function in this package; this runs them together, with the pieces a reviewer needs and that are easy to forget on a deadline: exact intervals on every change, multiple-comparison adjustment over the groups scanned, the envelope that rounding and suppression in the release imply, and a drift screen against the prior capsule.

Usage

anltab(
  data,
  value,
  prd,
  by = NULL,
  pop = NULL,
  prior = NULL,
  units = c("count", "continuous", "percent"),
  rnd = NULL,
  rndkind = c("nearest", "random"),
  suplim = NULL,
  dir = c("neutral", "higher", "lower"),
  per = 1000,
  conflvl = 0.95,
  alpha = 0.05,
  adjust = "BH",
  prdord = NULL
)

Arguments

Details

The steps, in order:

1. profile: profcols() and misssmry() on data;

2. change: yoy() between consecutive periods within each by group. For units = "count" the exact p-values are then adjusted across all comparisons with scanadj(), and, when rnd or suplim is given, yoybnd() adds the publication envelope (see pubbnd()). For "continuous" and "percent" units neither step runs: rnd and suplim are then ignored, with a warning;

3. rates (only with pop): rtcalc() per group and period, and ratechg() between periods, adjusted with scanadj();

4. trend: per group, the values summed within each period, then trndtst() (Mann-Kendall tau, p-value and Theil-Sen slope) and, for counts, cnttrnd() (Poisson rate ratio per period). A group with fewer than three periods gets NA; a test that fails gets NA rather than an error;

5. drift (only with prior): cpsdrft() of data against prior.

Nothing is written and nothing is fetched; a progress bar for the trend step is drawn on stderr unless options(morie.quiet = TRUE). Errors: data or prior not a data frame; a named column not in data; an invalid units, dir or rndkind; a conflvl outside (0, 1); and any error of the functions called (for example negative counts in yoy()).

Value

A list of class bricklayer_analysis (with a print method that summarises it in a few lines):

profile

profcols() of data.

missingness

named numeric from misssmry().

change

the yoy() table; for counts with the added columns p_value, padj and significant, and with env_* and compctlo/compctup when publication bounds were requested.

rates, ratechg

the rtcalc() and ratechg() tables when population is given, otherwise NULL.

trend

data frame with one row per group: the by columns, periods, first, last, tau, trndp, slope, rateratio, rateratiolo, rateratioup.

drift

cpsdrft() result against prior, or NULL.

meta

list of the arguments used, nrows, the sorted periods, prior (logical) and when (a time stamp).

See Also

rptanl() writes the result as a Markdown or HTML report.

Examples

path <- system.file("extdata", "otis_a01_individuals.csv",
                    package = "rmoriebricklayer")
otis <- read.csv(path)
a <- anltab(otis, value = "individuals", prd = "year",
                   by = c("table", "group"))
a
a$change[, c("table", "group", "year", "pctchg", "padj")]
a$trend[, c("group", "tau", "trndp", "rateratio")]

# counts published rounded to base 5: the intervals widen
b <- anltab(otis, value = "individuals", prd = "year",
                   by = c("table", "group"), rnd = 5)
b$change[1:3, c("group", "year", "pctlo", "compctlo")]

# rates with an exposure column
d <- data.frame(year = rep(2019:2022, 2), site = rep(c("A", "B"), each = 4),
                n = c(10, 12, 15, 9, 30, 28, 35, 40),
                pop = rep(c(1000, 2000), each = 4))
r <- anltab(d, "n", "year", by = "site", pop = "pop", per = 1000)
r$ratechg[, c("site", "year", "rate", "rateratio", "padj")]

# a column that is not in the data is an error
try(anltab(otis, value = "count", prd = "year"))

Description

Downloads each URL from here, records its SHA-256, size and time in a manifest, timestamps the manifest with an RFC 3161 time-stamping authority, and, with archive.org keys, uploads the files, the manifest and its timestamp as one public archive.org item. The item is the mirror that fetchvrf falls back to, and the timestamp shows that these were the bytes on that date.

Usage

arccap(
  urls,
  dir,
  item = NULL,
  meta = list(),
  tsa = "https://freetsa.org/tsr",
  ak = NULL,
  sk = NULL,
  tmout = 600L,
  upload = NULL,
  reason = NULL,
  visibility = c("public", "private"),
  warc = TRUE
)

Arguments

Details

Some servers refuse the Internet Archive's crawler while serving everyone else: Save Page Now then reports that the target “cannot find” the page, or “blocks access” to it, although the page is there. Those URLs cannot be put in the Wayback Machine, but they can be captured by anyone who can reach them, and the capture can be published. Pass a wbsave result as urls to capture just the URLs it failed on with such a refusal.

The manifest (manifest.json) lists, for each URL, the file, its SHA-256, its size and the time it was fetched. The timestamp (manifest.tsr) is the authority's signed statement that the manifest's SHA-256 existed at its time; tsvrf checks it, as does openssl ts -verify. Each file's digest in the manifest ties the file to that time.

The upload uses the Internet Archive's S3-like API: each file is sent with a PUT, streamed from disk, to
https://s3.us.archive.org/<item>/<file>.
The first request creates the item with the metadata given. The keys are the same as for wbsave, so they are read from IA_S3_ACCESS_KEY and IA_S3_SECRET_KEY, or their _FILE twins. They travel only in the authorization header and are never written or printed. Publish only what the source's licence lets you redistribute.

Value

A data frame of class rmbl_archive_capture, one row per URL: url, file, sha256, bytes, fetched_utc, status ("captured", "uploaded", "archived" or "failed: ..."), mirror (the archive.org download address, for fetchvrf), warc (the local WARC file), capture_id and archive (its page at archive.rmorie.com) and wayback (the Wayback Machine's answer); attributes manifest, timestamp (a tsvrf result, or NULL) and item.

References

Adams C, Cain P, Pinkas D, Zuccherato R (2001). Internet X.509 Public Key Infrastructure Time-Stamp Protocol (TSP). RFC 3161. DOI 10.17487/RFC3161

Internet Archive. ias3: the Internet Archive's S3-like API. archive.org/developers/ias3.html

See Also

fetchvrf, wbsave, tsvrf.

Examples

# save from your own connection to archive.rmorie.com and the Wayback Machine
# at once
src <- "https://www150.statcan.gc.ca/n1/tbl/csv/35100177-eng.zip"
cap <- arccap(
  src,
  file.path(tempdir(), "cap"),
  upload = c("rmorie", "wayback"),
  reason = "Custody table cited in chapter 3"
)
cap$archive
w <- wbsave(urls)
lic <- "https://www.statcan.gc.ca/en/terms-conditions/open-licence"
meta <- list(
  title = "Statistics Canada tables (captures)",
  creator = "Statistics Canada",
  licenseurl = lic
)
cap <- arccap(
  w,
  file.path(tempdir(), "capture"),
  item = "my-statcan-tables-2026",
  meta = meta
)
fetchvrf(
  cap$url[1],
  tempfile(),
  sha256 = cap$sha256[1],
  mirrors = cap$mirror[1]
)

Description

Lists the captures the archive holds of a web address, oldest first: when each was taken, its HTTP status, media type, size and SHA-256, and the addresses to open or download it. An address ending in * lists every address that starts with it. Public captures are listed for anyone; with your key (see dlogin) your private ones are listed too.

Usage

arcfind(url, from = NULL, to = NULL, limit = 1000L)

Arguments

Details

The archive answers the same CDX query as the Wayback Machine (/cdx/search/cdx?url=...&output=json), so tools written for it can read it as well. Set MORIE_ARCHIVE_URL to query another archive that speaks the same protocol, with its key in MORIE_ARCHIVE_KEY.

Value

A data frame of class rmbl_archive_find: capture_id, captured (UTC time), url, status, mime, bytes, sha256, mode ("client": saved from the capturer's connection; "worker": saved by the archive), visibility, capture (its page) and replay (its exact bytes).

See Also

arcget, arccap, fetchvrf.

Examples

caps <- arcfind("https://www150.statcan.gc.ca/n1/tbl/csv/35100177-eng.zip")
caps
arcfind("www150.statcan.gc.ca/n1/tbl/*")

Description

Downloads a capture's bytes and checks their SHA-256: against the archive's record, and against sha256 when you give one (the pin from your analysis). A file that does not match is deleted.

Usage

arcget(x, dest, sha256 = NULL, at = NULL, tmout = 600L)

Arguments

Value

dest, invisibly, with attributes capture_id, url, captured and sha256.

See Also

arcfind, fetchvrf (which falls back to the archive on its own).

Examples

src <- "https://www150.statcan.gc.ca/n1/tbl/csv/35100177-eng.zip"
f <- arcget(src, tempfile(fileext = ".zip"))
attr(f, "sha256")

Description

arccap() keeps every capture on this computer first: the files, their SHA-256, a manifest and, by default, a WARC copy. When the Wayback Machine, archive.rmorie.com or archive.org is down, or there is no network but the source's, the capture is still complete there. arcpush() sends that same capture on later: the very files that were fetched, each checked first against the SHA-256 the manifest recorded, so an archive keeps what you captured, not a new download. Run it again to retry what failed; what an archive already took is not sent twice.

The Wayback Machine takes no uploads: with upload = "wayback" it is asked to capture each address itself, as the page is then. If the page has changed since, its copy differs from yours, and your copy, with the manifest's timestamp, remains the record of what you had.

A manifest with no timestamp yet (the authority was unreachable, or tsa = NULL at capture) is stamped now when tsa is given; that timestamp shows the files existed by now, not by the time they were fetched.

Errors: dir without the ‘manifest.json’ that arccap() writes, or a manifest that lists no files; nowhere to send (neither item nor upload); as for arccap() for item, meta, upload, reason, tmout and the keys.

Usage

arcpush(
  dir,
  item = NULL,
  meta = list(),
  ak = NULL,
  sk = NULL,
  upload = NULL,
  reason = NULL,
  visibility = c("public", "private"),
  tsa = "https://freetsa.org/tsr",
  tmout = 600L
)

Arguments

Value

An object of class rmbl_archive_capture, as arccap() returns: one row per file of the capture, with what each archive answered, this time or an earlier one. A file that is missing, or whose bytes no longer match the manifest, is "failed: ..." and is not sent. What each archive took is also kept in ‘push.json’ in dir, which the next run reads.

See Also

arccap() makes the capture; wbsave() and wbrtry() for the Wayback Machine alone; fetchvrf() falls back to the copies.

Examples

# capture now, on this computer only
src <- "https://www150.statcan.gc.ca/n1/tbl/csv/35100177-eng.zip"
d <- file.path(tempdir(), "cap")
arccap(src, d, tsa = NULL)
# later, when the archives answer: the same bytes, sent on
arcpush(
  d,
  upload = c("rmorie", "wayback"),
  reason = "Custody table cited in chapter 3"
)

Description

Returns x unchanged when it is valid, well-formed text (so legitimate UTF-8 such as an accented name is preserved), and only transliterates to plain ASCII via toascii() when the text is not valid UTF-8 (an encoding error) or when force = TRUE (for ASCII-only destinations such as a package DESCRIPTION) . This lets author and supervisor names keep their accents wherever UTF-8 is supported while degrading gracefully instead of erroring where it is not.

Usage

asciifb(x, force = FALSE)

Arguments

Details

Each element is tested on its bytes with validUTF8(), so the outcome does not depend on the session locale. Valid elements (including plain ASCII and properly encoded accented text) are returned unchanged; invalid ones go through toascii(), which drops the bytes that are not UTF-8 and transliterates the rest. With force = TRUE every element goes through toascii(). A missing value is kept as NA unless force = TRUE, when it becomes "".

Value

A character vector the length of x: each element unchanged when it is valid UTF-8 (and force is not TRUE), otherwise its ASCII transliteration.

See Also

toascii(), which does the transliteration.

Examples

# By default valid UTF-8 is preserved (accents kept where supported).
asciifb("\u00c1ngela")               # "\u00c1ngela"

# force = TRUE always transliterates (for ASCII-only destinations
# such as a package DESCRIPTION).
asciifb("\u00c1ngela", force = TRUE)  # "Angela"

# Plain ASCII is returned unchanged either way.
asciifb("plain name")

# Vectorised; each element handled independently.
asciifb(c("caf\u00e9", "resume"), force = TRUE)

# Bytes that are not UTF-8 (a Latin-1 file read as UTF-8) are repaired.
bad <- rawToChar(as.raw(c(0x63, 0x61, 0x66, 0xe9)))
validUTF8(bad)
asciifb(bad)

Description

The share of all cases in a population that is due to an exposure: the fraction that would not occur if the exposure were removed and nothing else changed. It combines how strongly the exposure raises risk (the relative risk) with how common it is (the prevalence); multiplied by the number of cases it gives the attributable cases, as polbrdn does.

Usage

attrfrac(rr, expoprev, rrlo = NULL, rrup = NULL)

Arguments

Details

Levin's formula (Levin 1953; Rothman, Greenland and Lash 2008, chapter 5):

$PAF = p (RR - 1) / (1 + p (RR - 1)),$

with $p$ the proportion exposed. With everyone exposed it is $1 - 1/RR$. It assumes the relative risk is causal and unconfounded and that removing the exposure changes nothing else. A relative risk below 1 gives a negative fraction (the exposure protects).

When both confidence limits of the RR are given, the fraction's limits are the formula applied to them: the formula increases with the RR, so the interval of the fraction is the image of the interval of the risk. Refusals: an rr or limit that is not a single finite positive number (logical values included), a prevalence outside $[0, 1]$, only one of the two limits, or limits that do not enclose rr.

Value

A single number, the attributable fraction: in [0, 1) for rr >= 1, negative for a protective exposure. When the limits were given it carries an attribute ci, a named numeric vector (lower, upper). It has class c("rmbl_attrfrac", "rmbl_measure", "numeric"), so it prints as a box with its inputs (kept in attribute "lakner"); arithmetic, [, c() and as.numeric() give plain numbers.

References

Rothman, K. J., Greenland, S. and Lash, T. L. (2008). Modern Epidemiology, 3rd ed., chapter 5. Levin, M. L. (1953). The occurrence of lung cancer in man. Acta Unio Internationalis Contra Cancrum 9, 531-541.

See Also

crfpm25 for a relative risk from a concentration, polbrdn for the attributable cases.

Examples

attrfrac(1.5, 0.4)
attrfrac(2, 1)      # everyone exposed: 1 - 1/RR

# With the confidence limits of the relative risk.
paf <- attrfrac(1.08, 1, rrlo = 1.06, rrup = 1.09)
attr(paf, "ci")

# A protective exposure gives a negative fraction.
attrfrac(0.8, 0.5)

# A prevalence above 1 is refused.
try(attrfrac(1.5, 40))

Description

Checks every categorical column of a data frame for the coding problems that silently change results: numbers that are really category codes whose labels were lost, labels still in foreign value-label form, labels that differ only by case or by an invisible space, empty strings or "NA" stored as a category, unused levels, and identifiers mistaken for categories. Run it once after every import, before any model. It also reports which level R will use as the reference category.

Usage

audcats(data, cols = NULL)

Arguments

Details

For each column the levels are levels(x) for a factor and the sorted distinct values (radix sort, as factor() would order them) otherwise. A hazard is reported when:

• the column is haven labelled (it still carries value labels and must be decoded with deccodes or declbld);

• every level looks like a number (^[0-9.]+$);

• every level starts with a code prefix such as "1. White";

• a level has leading or trailing whitespace, including non-breaking and other Unicode spaces, or two levels are equal once trimmed;

• a level is empty, or is a missing-value sentinel (NA, N/A, NAN, NULL, NONE, ., -, ?, compared case-insensitively);

• the reference (first) level is empty, a sentinel, or padded;

• two different trimmed levels agree once lower-cased;

• a factor has levels that never occur;

• there are more than 50 levels.

The function only reports; it never changes data. It errors when data is not a data frame.

Value

A data frame of class c("bricklayer_category_audit", "data.frame"), one row per audited column, with character columns column, storage (the column's class), levels (up to eight levels joined with "|"), reference (the first level) and hazards (all hazards joined with " ;; ", or ""), and integer n_levels. The attribute clean is TRUE when no hazard was found. The print method lists each column and its hazards.

See Also

deccodes, grdlvls, xfrvrf.

Examples

# codes without labels, and a case-variant city name
a <- audcats(data.frame(
  race = factor(c("1", "2", "2", "3")),
  city = c("Toronto", "toronto", "Ottawa")[c(1, 2, 3, 3)]
))
a
attr(a, "clean")

# padded, empty and sentinel labels, and an unused level
audcats(data.frame(
  g = factor(c(" A", "A", "NA", ""), levels = c("", " A", "A", "NA", "B"))
))

# a clean frame; non-categorical columns are skipped by default
audcats(data.frame(sex = c("F", "M", "F"), age = c(31, 45, 27)))

try(audcats(1:3))

Description

Encodes and decodes base64, in both the standard alphabet and the URL-safe variant. Computed in R with no dependency, so a capsule can embed binary content in a text manifest wherever this package runs.

Usage

jsb64enc(input)

jsb64dec(input)

jsb64uenc(input)

jsb64udec(input)

Arguments

Details

jsb64enc() breaks its output into 72-character lines (a newline after each complete block of 54 input bytes), matching jsonlite's encoder byte for byte; the decoder ignores line breaks, other white space and = padding, so either form round trips. Any other character outside the alphabet is an error rather than being skipped, because skipping it would decode a corrupted body into a shorter valid-looking one.

The URL-safe variant substitutes - and ⁠_⁠ for + and / and drops the = padding, which is what makes it safe in a URL path, a query string or a filename.

Base64 is an ENCODING, not encryption or a digest: it hides nothing and anyone can reverse it. Use csha256() to pin content and chmac() to authenticate it.

Errors: an encoder given anything but a raw vector, a character vector or NULL; jsb64dec() given anything but text or raw bytes, or text with a character outside the standard alphabet; jsb64udec() given anything but character text (raw bytes are refused), NA, or a character outside the URL-safe alphabet.

Value

The encoders return a character scalar ("" for empty input, NA_character_ for NULL); the decoders return a raw vector (use rawToChar() for text).

See Also

jsgzenc(), which composes this with gzip.

Examples

# Round trip through the standard alphabet.
b <- jsb64enc("hello capsule")
b
rawToChar(jsb64dec(b))

# Raw input works the same way.
jsb64enc(charToRaw("abc"))
jsb64dec(jsb64enc(charToRaw("abc")))

# Padding appears when the length is not a multiple of three.
jsb64enc("a")
jsb64enc("ab")
jsb64enc("abc")

# The URL-safe variant has no "+", "/" or "=" to escape.
jsb64uenc(as.raw(c(255, 224, 63)))
jsb64enc(as.raw(c(255, 224, 63)))
jsb64udec(
  jsb64uenc("path/safe?yes")
)

# Long input is wrapped, and the decoder ignores the breaks.
long <- jsb64enc(strrep("x", 200))
grepl("\n", long)
rawToChar(jsb64dec(long)) == strrep("x", 200)

# A character outside the alphabet is refused, not skipped.
try(jsb64dec("ab$c"))

Description

Splits the two-way fixed-effects difference-in-differences estimate of a staggered adoption into the two-group, two-period comparisons it averages, with their weights: treated against never treated, earlier against later treated, and later against earlier treated (the comparisons that use already-treated units as controls).

Usage

bacondcmp(data, unit, time, y, trtd)

Arguments

Details

For a balanced panel with treatment that stays on and no covariates, the estimate of plndid is exactly $\hat\beta = \sum_{k} s_{kU}\hat\beta_{kU} + \sum_{k \ne l} s_{kl}\hat\beta_{kl}$, with weights that sum to one and grow with the groups' sizes $n$ and the variance of their treatment over time $\bar D(1 - \bar D)$ (Goodman-Bacon 2021, theorem 1). Units treated from the first period act as an already-treated group. A large weight on “later vs earlier” comparisons, together with effects that grow over time, is the warning sign that the two-way estimate is biased.

Value

A data frame of class rmbl_bacon: type ("treated vs never", "earlier vs later", "later vs earlier"), treated, control (the adoption periods; Inf for never), estimate, weight; attribute twfe (the weighted sum, equal to the two-way estimate).

References

Goodman-Bacon A (2021). Difference-in-differences with variation in treatment timing. Journal of Econometrics 225(2):254-277. DOI 10.1016/j.jeconom.2021.03.014

See Also

plndid, stgatt.

Examples

p <- polex()
b <- bacondcmp(p, "district", "month", "crmrate", "reformed")
b
tapply(b$weight, b$type, sum)

# your own panel: eight towns over ten years, four adopting a by-law (two
# in 2019, two in 2021), with each town's population, county and income
own <- data.frame(
  town = rep(c("A", "B", "C", "D", "E", "F", "G", "H"), each = 10),
  year = rep(2015:2024, times = 8)
)
i <- match(own$town, unique(own$town))
own$first <- c(2019, 2019, 2021, 2021, NA, NA, NA, NA)[i]
own$bylaw <- as.integer(!is.na(own$first) & own$year >= own$first)
own$pop <- c(12, 30, 8, 22, 15, 40, 9, 18)[i] * 1000
own$county <- rep(c("North", "South"), 4)[i]
own$income <- 50 + own$year - 2015 + cos(seq_len(80))
own$rate <- 30 + i + 0.5 * (own$year - 2015) - 4 * own$bylaw +
  0.2 * own$income + sin(seq_len(80))
bacondcmp(own, "town", "year", "rate", "bylaw")

Description

Means of each covariate in the treated and comparison groups, their standardised mean difference, and Welch's test, optionally weighted (after matching or weighting).

Usage

baltab(data, grp, covs, w = NULL)

Arguments

Details

The standardised mean difference is $(\bar x_1 - \bar x_0) / s_1$ with $s_1$ the unweighted standard deviation in the treated group (the convention for the effect on the treated, so that it does not move when weights change the comparison group); $|SMD| < 0.1$ is the usual threshold for balance. With weights the means are weighted and the test uses the weighted variances and Kish's effective sample sizes.

Value

A data frame of class rmbl_balance: covariate, meantrtd, meanctl, smd, p_value.

References

Welch BL (1947). The generalization of ‘Student's’ problem when several different population variances are involved. Biometrika 34(1-2):28-35. DOI 10.1093/biomet/34.1-2.28

Examples

d <- data.frame(t = rep(0:1, each = 50),
                age = c(stats::rnorm(50, 30), stats::rnorm(50, 33)))
baltab(d, "t", "age")

Description

Recomputes a statistic of banded data (by default the Gini coefficient) for a range of assumed caps on the open top band, and reports how far the answer moves. Every statistic computed from banded data depends on where the open band is assumed to end; this function measures that dependence, so you can report a range instead of a single figure when the dependence is large.

Usage

bandsens(bands, cnts, stat = ginicoef, caps = NULL, rule = "midpoint")

## S3 method for class 'rmbl_band_sensitivity'
format(x, digits = 4L, ...)

## S3 method for class 'rmbl_band_sensitivity'
print(x, ...)

Arguments

Details

Let $b$ be the smallest lower bound among the open top bands. The default caps are 12 points spaced evenly on the log scale from $1.05 b$ to $20 b$. For each cap $c$, bandvals(bands, rule, openupcap = c) gives one value per band, each value is repeated cnts times, and stat is applied to the result. A statistic that errors or does not return a single number gives NA; warnings are suppressed.

The span is $\max - \min$ over the finite values and the relative span is the span divided by the median absolute value (NA when any value is zero). The print method adds a note when the relative span exceeds 0.1.

Errors: cnts not numeric or of a different length from the bands; no open top band.

Value

A data frame of class c("rmbl_band_sensitivity", "data.frame") with numeric columns cap and value, and attributes span, relspan and base (the lower bound $b$ of the open band).

Examples

b <- prsbands(c("1", "2 to 5", "6 to 10", "Greater than 10"))
counts <- c(1200, 430, 110, 38)

s <- bandsens(b, counts)
s
attr(s, "relspan")

# A statistic that barely moves has been measured; one that swings
# has not.
bandsens(b, counts, stat = mean, caps = c(12, 20, 50))

# nothing to vary without an open band
try(bandsens(c("1", "2 to 5"), c(10, 4)))

Description

Gives each band a single number that a calculation can use (a mean, a concentration index, a regression), under a rule you choose. For closed bands the rule places the value at the midpoint, the lower or upper bound, or the geometric mean of the bounds. An open band ("50+", “under 18”) has no midpoint, so its missing bound is filled from an explicit assumption, and the rows that rest on that assumption are flagged in the result instead of being absorbed into it.

Usage

bandvals(
  bands,
  rule = c("midpoint", "lower", "upper", "geometric"),
  openupcap = NULL,
  openupfac = 2,
  openloflr = NULL
)

Arguments

Details

For a band with bounds $[l, u]$ the value is $(l + u)/2$ ("midpoint"), $l$ ("lower"), $u$ ("upper") or $\sqrt{l u}$ ("geometric"; when $l \le 0$ the arithmetic midpoint is used instead, because the geometric mean is undefined there). The "lower" and "upper" rules bracket any value that could come from the bands.

An open upper band gets $u =$ openupcap if given, else $u = l \times$ openupfac. An open lower band gets $l =$ openloflr if given, else $\min(0, u)$. A cap below the band's lower bound, or a floor above its upper bound, is an error, because the value would lie outside the band. Unparsed bands give NA.

Use bandsens to see how much a result moves with the cap, and expbands to turn the values into one value per unit.

Value

The band table with two added columns: value (numeric, the representative value) and assumed (logical, TRUE for the open bands whose value rests on the cap or floor).

See Also

prsbands, bandsens, expbands

Examples

b <- prsbands(c("1", "2 to 5", "6 to 10", "Greater than 10"))
bandvals(b)

# The lower and upper rules bracket the truth.
bandvals(b, rule = "lower")$value
bandvals(b, rule = "upper", openupcap = 40)$value

# geometric placement, and an explicit floor for an open bottom band
bandvals(b, rule = "geometric", openupcap = 40)$value
bandvals(c("under 5", "5 to 9", "10+"), openloflr = 1)

# a cap inside the open band is refused
try(bandvals(b, openupcap = 5))

Description

Compares the distribution of leading significant digits in x with Benford's law, $P(d) = \log_{10}(1 + 1/d)$. Naturally occurring quantities that span several orders of magnitude follow it closely; figures that were rounded, truncated, capped, re-scaled, or invented typically do not. That makes it a cheap screen for a numeric column that arrived looking plausible but is not the measurement it claims to be.

It is a SCREEN, not a verdict. Columns with a narrow range, a unit floor or ceiling, or an assigned-identifier structure (postcodes, year fields, prices ending in 99) legitimately violate Benford's law. Treat a small p-value as a reason to look, never as evidence of fabrication.

Usage

benfordtst(x)

Arguments

Details

Zeros and non-finite values have no leading significant digit and are excluded; the sign is ignored.

With $n$ usable values, $O_d$ the count with leading digit $d$ and $E_d = n \log_{10}(1 + 1/d)$,

$X^2 = \sum_{d=1}^{9} \frac{(O_d - E_d)^2}{E_d}$

on 8 degrees of freedom. The smallest Benford share is 0.046 (digit 9), so the chi-square approximation needs about $n \ge 110$ for every expected count to reach 5; below that treat the p-value as rough. Digits are counted in compiled code in blocks of a million values, and a long vector shows a progress bar in an interactive session (silenced by options(morie.quiet = TRUE)). Errors: non-numeric input, or no finite non-zero value.

Value

A list of class bricklayer_benford with counts (observed digit frequencies 1–9), expected, proportion, statistic, df, p_value, logpvalue (the log of the upper tail, which keeps its information where p_value underflows to 0), and n (the number of values used). counts, expected and proportion are named numeric vectors with names "1" to "9". The print method shows the verdict and a bar chart of observed against expected shares.

References

Benford F (1938). The law of anomalous numbers. Proceedings of the American Philosophical Society 78(4), 551–572.

See Also

drftchisq(), cpsdrft().

Examples

# A quantity spanning several orders of magnitude follows the law.
set.seed(3)
benfordtst(10^stats::runif(2000, 0, 6))

# Digits drawn uniformly do not.
benfordtst(as.numeric(paste0(sample(1:9, 2000, TRUE), "000")))

# The expected proportions are the closed form.
b <- benfordtst(10^stats::runif(500, 0, 5))
all.equal(unname(b$expected) / b$n, log10(1 + 1 / (1:9)))

# Zeros carry no leading digit and are excluded from n.
benfordtst(c(0, 0, 1, 2, 3))$n

try(benfordtst(c(0, NA)))

Description

Dispatches the verbs of the rmoriebricklayer launcher:

login [--token] [--email ADDRESS [--code CODE]] [--no-browser]

sign in for the curated tables at data.rmorie.com (a free MORIE data key, issued at www.rmorie.com/access/): a key you paste with --token (read from the terminal or a pipe when KEY is omitted), a code sent to --email, or the GitHub device flow

logout

forget the data key

doctor

report the language-model routes available here (own endpoint, each provider with your key, local Ollama), in the order ask tries them

models [--route ROUTE]

list the models your provider (or server) offers, the one in use marked

data list

the curated tables at data.rmorie.com

data pull db/table

download one table as CSV

app [PATH] [--port N] [--no-browser]

open the analysis dashboard (rmblapp) on this machine, with PATH, a CSV file or a capsule directory, loaded; needs the shiny package

config [setup | set KEY VALUE | unset KEY | get KEY | help]

show or change the settings of the language model (llmcfg): the route ask uses and the address, key and model of each route; setup walks through them

ask [--route ROUTE] [--model NAME] [--fresh] PROMPT...

send a prompt to the model (or the named one, on the named route) and print the reply; the same question again is answered from this machine (--fresh asks anew). aye and tellme are the same verb, and a quoted question with no verb (rmbl "how do I weight a survey?") is asked too

explain-error MESSAGE...

where an error of these packages comes from (explerr), without a model or the network

bundle REQUEST...

agentbndl from the shell

pack DIR [--dest D] [--data FILE] [--version V] [--overwrite] [--secrets]

a reproducibility bundle of the project in DIR, or of a shipped example such as otis-mrp, with mkbndl: one zip a reviewer unzips and double-clicks (default --dest: the working directory; --secrets packs files that look like credentials, left out otherwise); pack --new DIR [--title T] [--author A] [--contact C] writes a starter project with mkbndlproj; pack --check ZIP [--data FILE] [--synthetic] [--keep] [--container[=podman|docker]] runs a built bundle as a recipient would and checks its results with mkbndlchk (exit status 1 when it fails). It RUNS THE BUNDLE'S CODE on this machine, with a minimal environment and a temporary HOME: check only bundles you trust. --container runs it in a podman or docker container; otherwise a PID namespace (Linux), a job object (Windows) or only a process group (macOS) contains it

wayback FILE [--out FILE.jsonl] [--jobs N] [--within 30d] [--max N]

save every URL listed in FILE (one per line; blank lines and # comments skipped) in the Internet Archive's Wayback Machine with wbsave, results in --out (default FILE.wayback.jsonl), so a stopped run carries on; wayback --retry FILE.jsonl [--include-gone] retries the failures of an earlier run (wbrtry)

archive save URL... --reason TEXT [--dir D] [--private] [--also wayback]

save from this connection to archive.rmorie.com with arccap, keeping the files, their manifest and WARC copies in D (default ‘captures’); archive push DIR [--reason TEXT] [--private] [--also wayback] sends a capture kept in DIR on later, the same bytes, with arcpush (run it again to retry what failed); archive find URL and archive get ID|URL -o FILE list the captures of a URL and fetch one's exact bytes (arcfind, arcget)

verify PATH [--json] [--fail-on LEVEL] [--public-key HEX]

check every capsule under PATH as the rmbl GitHub Action does (cpschk): each file against its JSON Schema (schemapath), vrfcps and cpsbndlvrf. The secret of an HMAC bundle is read from RMBL_KEY, which is then removed from the environment. Prints the Markdown report, or JSON with --json; exits 1 when a failed row is at or above --fail-on (error, the default, drift or warn). --base MANIFEST|DIR compares each manifest with its base copy, --manifests GLOB checks only the capsule files matching GLOB, --snapshot looks up each source's latest Wayback snapshot and --no-footprint leaves the footprint (co2fp) out of the report. --seen-log FILE adds the signatures an XMSS seen-log lists to those under PATH: an XMSS index one key signed twice fails its xmss_index row

drift PATH [--json] [--fail-on LEVEL] [--html FILE] [--app]

verify, then fetch each source again and compare it column by column (cpsdrft, Benford): a column past --psi-warn X (0.1) is a warning, a flagged one past --psi-fail X (0.25) drift. --allow-file accepts file:// sources, --cache DIR keeps each source and reads it when a fetch fails; the options of verify apply. --html FILE writes the interactive report (drfthtml) and --app opens it as a Shiny app (drftapp)

refresh PATH

pin each source's new bytes (cpsrfrsh) and print what changed (mnfdiff)

report PATH [--format md|json]

the verify report, exit 0; the options of verify apply

sign PATH [--scheme S] [--key-file F] [--overwrite]

pin each provenance on its first run and sign each capsule into capsule_bundle.json (cpsseal), then verify as verify does. The key is RMBL_KEY, removed from the environment at once, or the contents of --key-file, and is never printed: the shared secret for hmac (the default), else the JSON of a fipskg key, or its secret hex with --public-key, or the JSON of a pqckg key for xmss-sha256 (see the section XMSS keys). An existing bundle is refused without --overwrite

functions [PATTERN]

list the exported functions with their one-line titles, optionally filtered by a regular expression

describe NAME

title, usage, arguments and description of one function, from its help page

examples NAME

print the examples of one function

version

print the package version (also --version, -v)

help [VERB]

this list, or one verb's usage (also --help, -h; VERB --help shows a verb's usage without running it)

The launcher installed by instcli calls this function with the shell's arguments and exits with its status. Calling it from R is the same as typing the command, which is how the examples below work.

Usage

blcli(args = commandArgs(trailingOnly = TRUE), out = cat)

Arguments

Details

An option a verb does not take is refused with exit status 2 rather than ignored or sent as part of a prompt; everything after a bare -- is positional (for a prompt that itself starts with --) and only ask, bundle and explain-error accept it. A first argument that is plainly a question (two or more words, or ending in ?, with no path, line break or key-shaped word) and no verb is asked as ask would ask it.

Errors are not raised: every error is written to out as "PROG: message" and turned into the exit status. The status is 0 on success, 1 when the verb failed (an error, nothing found, no route), 2 for a usage error (an unknown verb or option, a missing value) and 130 when interrupted (Ctrl-C). The command name in messages is the environment variable RMBL_PROG, else rmoriebricklayer.

Side effects: during the call options(morie.progress = TRUE) is set, so downloads draw a progress bar; when RMBL_PROG is set (the launcher sets it) and options(rmoriebricklayer.dcache) is not, the data.rmorie.com manifest and tables are kept in tools::R_user_dir("rmoriebricklayer", "cache") rather than a fresh temporary folder. The verbs login, logout and config write the shared credentials and settings files; ask, models, data, bundle and wayback use the network, and drift and refresh fetch each capsule's source (refresh writes the capsule), sign writes each capsule's provenance and bundle; verify and report only read the capsule files; functions, describe, examples, explain-error, version and help only read the installed help.

Value

The exit status (integer scalar), invisibly: 0 on success, 1 on failure, 2 for a usage error, 130 when interrupted.

XMSS keys

An XMSS key (pqckg) signs each of its 2^height indices once; two signatures at one index let anyone forge a third. Its private key never changes, but its state, the indices already used, changes with every signature, and a copy of an old state signs used indices again (RFC 10033, “Hash-Based Signatures: State and Backup Management”, sections 2 and 4). rmbl sign --scheme xmss-sha256 --key-file KEY therefore treats the key file as the state:

• Before anything is signed it rewrites KEY, under a lock, through a temporary file in the same directory that is flushed to the disk and renamed over it, with next_index advanced and the state in a used field: highest (the highest index used), generation (one more at every write), writer (a digest of the host name and boot id), the last reserved block and, once split, sectors. A crash after that write loses the reserved indices and never reuses one. A KEY that cannot be rewritten (not writable, a read-only file system) is refused.

• --state-file FILE keeps the state in FILE (an absolute path to a writable regular file, or a new file in a writable directory) and only reads KEY; a key in RMBL_KEY needs it. The command line never uses the configuration directory that cpssign uses in an R session; it reads and updates the record in RMBL_XMSS_INDEX_FILE only when that is set.

• --sectors K --sector NAME (K a power of two, at most half the indices) splits the indices into K disjoint ranges and signs only in one of them, which keeps its own highest index; the range of NAME is SHA-256 of NAME and the key fingerprint, modulo K. Two names land in one sector by chance, so for real use give each runner its own --sector-index I. A sector that is used up is refused. Once split, sign only through rmbl sign.

• --reserve N marks the next N indices used in one write (at least one per capsule signed); the capsules of the run are signed from that block, and the indices it does not use are lost, as RFC 10033, section 5.9 allows.

• A stale copy is refused: when a bundle under PATH, or the --seen-log FILE that every rmbl sign appends its signatures to, shows this key signed an index above what the state says, nothing is signed. A note says when another host wrote the state last.

• When N or fewer indices are left in the range (--warn-left N, by default a tenth of the range) a warning is printed in the report and on standard error (RFC 10033, section 3.4).

A CI secret is a fixed copy of the key that every job starts from, so each job would sign the same indices again: in CI sign with ML-DSA (fipskg), or keep the XMSS state in a file that persists between jobs and give each concurrent runner its own sector. rmbl verify shows each XMSS signature's index and key fingerprint, and xmssidxseen lists indices signed twice. NIST SP 800-208 approves XMSS only when the key is made and used inside a hardware cryptographic module that never exports it; a key file is software, so this follows the software guidance of RFC 10033 and is not a validated SP 800-208 deployment.

See Also

instcli to put the command on your PATH.

Examples

blcli("version")
blcli("help")

# one verb's usage
blcli(c("help", "data"))

# capture the output instead of printing it
con <- textConnection("lines", "w")
status <- blcli(c("functions", "^drbg_"), out = function(x) {
  cat(x, file = con)
})
close(con)
status
lines

# an unknown verb is a usage error: status 2, nothing run
blcli("nonsense")

# check a capsule as the rmbl GitHub Action does
d <- usecpstmpl(tempfile("capsule-"), ex = TRUE)
invisible(cpsrfrsh(d))
blcli(c("verify", d))


# these read every help page, which takes a few seconds
blcli(c("functions", "json"))
blcli(c("describe", "anltab"))
blcli(c("examples", "perdays"))

Description

Downloads one file from the web to a local path, showing a progress bar (percent, size and rate) while a person is watching. It is the single download routine of the package and its sibling packages, and it applies the package's network safety rules: HTTPS only, no local or private addresses, a size cap, and checked redirects. Use it when a script must fetch a source file and you want a visible, safe transfer.

Usage

bldl(
  url,
  dest,
  headers = NULL,
  label = basename(dest),
  size = NULL,
  tmout = 3600,
  quiet = NULL,
  tty = NULL,
  allowfile = FALSE,
  maxbytes = 2^31
)

Arguments

Details

Progress is an animated bar, redrawn in place, when the console can redraw a line: a terminal, or the RStudio, Positron, R.app or RKWard console (never a knitr or Jupyter notebook; the options rmoriebricklayer.dynamic or cli.dynamic, or the environment variable R_CLI_DYNAMIC, override the choice). In the RStudio or Positron console the bar is drawn on stdout, elsewhere on stderr. Otherwise milestone lines go to stderr (every 10 percent, or every 50 MB when the size is unknown), and nothing is shown when quiet. By default output is shown in an interactive session or when options(morie.progress = TRUE) is set, and suppressed by options(morie.quiet = TRUE).

Network: the transfer runs in the package's compiled HTTP client. The URL is checked, its host resolved and pinned, and every redirect is checked again (no switch from https to http; credentials are not passed to another host). Plain http:// is allowed only when rmoriebricklayer.allowhttp, an option, is TRUE. A transfer slower than 64 bytes a second for 30 seconds is abandoned. A file:// URL is copied locally through the same progress display, only with allowfile = TRUE.

Errors: a URL that is not https (or allowed http/file), or that points at a loopback, link-local or private address; headers not a named character vector; tmout below 1 or above 2147483647; maxbytes below 1; a network failure; an HTTP status of 400 or more; a body larger than maxbytes. dest is overwritten. A server busy for a moment (HTTP 429, 500, 502, 503 or 504) is asked again up to three times, after the wait its Retry-After header asks for (at most a minute) or else 1, 2 and 4 seconds.

Value

dest, invisibly.

See Also

blfetch for a download with a Wayback Machine fallback.

Examples

# offline: a local file through the same routine
src <- tempfile(fileext = ".txt")
writeLines("hello", src)
dest <- tempfile(fileext = ".txt")
bldl(paste0("file://", src), dest, quiet = TRUE,
                    allowfile = TRUE)
readLines(dest)

# the safety checks run before any connection is made
try(bldl("http://example.org/table.csv", tempfile()))
try(bldl("https://127.0.0.1/table.csv", tempfile()))
try(bldl(paste0("file://", src), tempfile()))

# a real download needs the network


bldl("https://www.ontario.ca/some/open-data/table.csv",
                    "table.csv", maxbytes = 50e6)

Description

The shared data-fetch foundation of the morie ecosystem. Downloads url to dest; if the live download fails (404, network), it retries the Wayback Machine snapshot – so a rotated or removed source file (CIHI, open-data portals) stays retrievable. Backed by the package's C++ kernel rmbl_fetch_with_fallback (libcurl), which rmorie and morie share through LinkingTo. The return value says which source served the file, so a report can note that it used an archived copy.

Usage

blfetch(url, dest, wb = "", tmout = 120L)

Arguments

Details

Before any network use the arguments are checked: url must be one https URL to a public host (plain http only when the option rmoriebricklayer.allowhttp is TRUE; loopback, private, link-local and metadata addresses are refused), the directory of dest must exist and be writable, and tmout must be at least one second. The compiled layer then resolves the host, pins the connection to the checked addresses and re-checks every redirect.

If the live download fails, the fallback is the wb URL when given, else a snapshot looked up through the Wayback availability API (as wbsnpnat() does). dest is overwritten. A spinner is drawn on stderr while waiting, unless options(morie.quiet = TRUE).

Errors: invalid arguments as above; and, when both sources fail, an error saying whether the live URL could not be reached or answered an HTTP status, and whether the Wayback Machine had no snapshot or its copy could not be downloaded.

Value

Invisibly, "live" or "wayback": the source that served the file. Throws on total failure.

Examples

# Inputs are validated before any network access:
try(blfetch("", tempfile()))          # empty url -> error
try(blfetch("ftp://example.org/data.csv", tempfile()))
try(blfetch("https://127.0.0.1/data.csv", tempfile()))
try(blfetch("https://example.org/data.csv",
                     file.path(tempdir(), "no-such-dir", "data.csv")))


# needs the network
dst <- tempfile(fileext = ".xlsx")

# Live download; auto-resolves a Wayback snapshot only if the live URL fails.
blfetch(
  paste0(
    "https://www.cihi.ca/sites/default/files/document/",
    "hospital-beds-2024-2025-data-tables-en.xlsx"
  ),
  dst
)

# Pin an explicit Wayback snapshot to fall back to, and a shorter timeout.
blfetch(
  "https://example.org/rotated-file.csv",
  tempfile(fileext = ".csv"),
  wb = paste0(
    "https://web.archive.org/web/2024id_/https://example.org/rotated-fil",
    "e.csv"
  ),
  tmout = 60
)

# The return value tells you which source served the file.
status <- blfetch("https://cloud.r-project.org/", tempfile())
status # "live" or "wayback"

Description

Every plot() method in the package draws with base graphics and the same small set of style choices, which these functions expose. blpal() returns colour-blind-safe palettes: "okabe" (the Okabe-Ito categorical set), "diverging" (blue, grey, orange), "sequential" (light to dark blue) and "grey". plotopt() sets the defaults every method uses, for the session: palette, grid, cex, las and bg. Each method also takes main, col (which overrides the palette), palette and ... passed to the underlying base-graphics call, so a single plot can be changed without changing the session. Use blpal() also to colour your own plots so they match the package's.

Usage

blpal(name = c("okabe", "diverging", "sequential", "grey"), n = NULL)

plotopt(...)

Arguments

Details

The palettes: "okabe" is the eight-colour set of Okabe and Ito (2008), distinguishable under the common forms of colour-vision deficiency. "diverging" has three anchors (blue #2166AC, grey #8C8C8C, orange-red #D6604D) and is used for changes (decrease, no change, increase). The palette "sequential" has five blues from light to dark; "grey" five greys from dark to light. With n given, "okabe" and "grey" repeat their colours (rep_len()) and "diverging" and "sequential" are interpolated with grDevices::colorRampPalette().

The session options are stored in options(rmoriebricklayer.plot); the defaults are palette = "okabe", grid = TRUE, cex = 1, las = 1, bg = "white" and panel = "#EEF1F6". Setting some options keeps the others.

Text is drawn at the session cex. Titles, axis titles, row labels and keys are made smaller to fit a small device, but never below 0.6 of that size, so that they stay readable; for larger text set, for example, plotopt(cex = 1.2), and for a table with many rows use a taller device.

Errors: an unknown palette name; n not a positive whole number; an option name other than the six above, or an unnamed value; an unknown palette option.

Value

blpal(): a character vector of hexadecimal colours. plotopt(): a named list of the six options in force (after the change), invisibly when setting.

References

Okabe, M. and Ito, K. (2008). Color Universal Design (CUD): how to make figures and presentations that are friendly to colorblind people. jfly.uni-koeln.de/color

Examples

blpal("okabe", 3)
# plotopt() returns the options in force, so keep them before a change
old <- plotopt()
plotopt(palette = "grey", grid = FALSE)
plotopt(palette = old$palette, grid = old$grid)

# an interpolated sequential scale, shown as swatches
cols <- blpal("sequential", 7)
barplot(
  rep(1, 7),
  col = cols,
  border = NA,
  axes = FALSE,
  names.arg = cols,
  las = 2,
  cex.names = 0.6
)

# the same palette in a base plot of your own
plot(1:8, pch = 19, cex = 3, col = blpal("okabe"))

# names are checked
try(blpal("viridis"))
try(plotopt(colour = "red"))

Description

Makes a progress bar for a loop or a long computation: the animated bar of bldl(), redrawn in place at no more than 12 frames a second, when the console can redraw a line – a terminal, or the RStudio, Positron, R.app or RKWard console (never a knitr or Jupyter notebook; the options rmoriebricklayer.dynamic or cli.dynamic, or the environment variable R_CLI_DYNAMIC, override the choice). In the RStudio or Positron console the bar is drawn on stdout, elsewhere on stderr. Otherwise a milestone line goes to stderr at each tenth, and nothing is shown when quiet. With total = 0 the work has no count, and the bar is a spinner with the time so far and the current stage (the info text), with a line at each new stage when the console cannot redraw.

Usage

blprog(label, total, quiet = NULL, color = "blue", delay = NULL)

Arguments

Details

Progress is shown when a person is watching: an interactive session, or a command line that sets options(morie.progress = TRUE). options(morie.quiet = TRUE) or quiet = TRUE silences it.

Nothing is drawn until the work has run for delay seconds, so a fast call never flashes a bar. A bar made inside a function while a bar made by a calling function is still running stays silent, so nested work shows one bar, for the outer loop.

Updates are cheap: when quiet an update only stores the value, and when drawing, frames are throttled, so $tick() can be called on every iteration of a loop. Call $done() through on.exit(p$done(), add = TRUE) right after making the bar, so the line is cleared and the bar is released even when the work fails or is interrupted:

p <- blprog("bootstrap", B)
on.exit(p$done(), add = TRUE)
for (b in seq_len(B)) {
  ...
  p$tick(info = sprintf("resample 
}

$done() may be called more than once; only the first call acts.

Value

An object of class rmbl_progress, a list of functions: update(done, info = NULL) sets the steps done (and the text after the bar), tick(by = 1, info = NULL) adds by steps, wait(seconds) sleeps while keeping the animation moving, and done(summary = NULL) ends the bar, clears the line and, when the bar was shown, writes summary on its own line. Each returns the object invisibly.

Examples

p <- blprog("squares", 50, delay = 0)
for (i in 1:50) {
  x <- sum(seq_len(i)^2)
  p$tick(info = sprintf("step %d/50", i))
}
p$done("squares: 50 steps")

# work with no count: a spinner and the stage
f <- function() {
  p <- blprog("model", 0, delay = 0)
  on.exit(p$done(), add = TRUE)
  p$update(info = "fitting")
  p$wait(0.2)
  p$update(info = "writing")
  invisible(TRUE)
}
f()

Description

Assigns each place (a block, an address) to the treated group when it lies within a radius of a site (a store that donated, a new facility), to the comparison group when it lies within the same radius of a comparison point and farther than a buffer from every site, and excludes it otherwise, so that spillover from the sites does not reach the comparison.

Usage

bufasgn(x, y, sitex, sitey, comx, comy, radius, buffer = 3 * radius)

Arguments

Details

Distances are Euclidean, so coordinates must be projected (metres or feet, not degrees). The radius and the buffer are in the same unit. The buffer should be several times the radius: the design it follows used 330 feet and 990 feet.

Value

A data frame of class rmbl_buffer with the columns group ("treated", "comparison" or "excluded"), site_distance and comparison_distance.

References

Applied to blocks near donor stores: Zorro Medina, Hackett, Green and Vargas (2025), see evstudy.

Examples

set.seed(4)
x <- runif(200, 0, 5000)
y <- runif(200, 0, 5000)
a <- bufasgn(x, y, c(1000, 3500), c(1200, 3000), c(2500, 4500), c(500, 4000),
                   radius = 330)
table(a$group)

Description

Records, for every package named, its version and enough about its provenance to find the same one again: the library it was loaded from, the repository it was installed from, and – for a package installed from a remote – the remote's URL and commit. Store it in a manifest next to capenv() when an analysis depends on packages installed from GitHub or a local build.

Usage

capdeps(pkgs = loadedNamespaces())

Arguments

Details

Why versions alone are not enough. Two installations can report the same version and differ: one built from CRAN, one from a fork, one from a local R CMD INSTALL of a working tree with uncommitted changes. A version number identifies an intention; the repository and commit identify what was actually loaded.

The fields are read with utils::packageDescription() from each package's installed DESCRIPTION: Repository (for example "CRAN"), RemoteUrl and RemoteSha (written by remotes, pak and devtools for packages installed from a Git host) and Built. Packages are sorted and de-duplicated. A progress bar is drawn on stderr for many packages unless options(morie.quiet = TRUE). Nothing is written. Error: pkgs that is not a character vector.

Value

A data frame with one row per package and character columns package, version, library (the library directory), repository, remote_url, remote_sha and built (R version, platform and time of the build). Unavailable fields are NA rather than omitted, so a reader can see that the information was absent rather than forgotten; for a package that is not installed, every field except package is NA, and a warning names it. An empty pkgs gives a data frame with no rows and the same columns.

See Also

capenv(), envdiff().

Examples

deps <- capdeps(c("stats", "utils"))
deps[, c("package", "version")]

# a package with no repository field records that fact
is.na(capdeps("stats")$repository)

# a package that is not installed is a row of NAs, not an error
capdeps("notapackage123")[, c("package", "version")]

# names must be character
try(capdeps(1))

Description

Records the facts a replicator needs to rebuild the session in which a result was produced: the R version, platform and operating system, a UTC timestamp, the random number generator in use and whether it was seeded, and the versions of the packages involved. Store the result in a results manifest (for example with wrmnfjs()) next to the numbers it produced, and compare two captures with envdiff().

Usage

capenv(pkgs = loadedNamespaces())

Arguments

Details

Package versions come from utils::packageVersion(); a package that is not installed gets NA rather than an error. The packages are sorted and de-duplicated. The generator kinds come from RNGkind() (R has changed its default sample() algorithm before, so a recorded seed means nothing without the kind it was fed to), and rng_seeded says whether .Random.seed exists in the global environment. Nothing is written and the random state is not touched.

Value

A list of seven elements:

r_version

character, e.g. "4.5.1".

platform

character, R.version$platform.

os

character, the system name and release.

captured_utc

character, ISO 8601 time YYYY-MM-DDThh:mm:ssZ.

rng_kind

character, the three RNGkind() values joined by commas.

rng_seeded

logical, whether .Random.seed exists.

packages

named list of version strings, sorted by name (NA for a package that is not installed). A list, not a vector, so the names survive the manifest's JSON: a named vector is written as a bare array.

The list has class c("rmbl_capenv", "list"), which only changes how it prints (a box); its elements, and the JSON written from it, are those of a plain list.

Examples

# Record specific packages' versions alongside the session facts.
env <- capenv(c("stats", "utils"))
env$r_version
env$os
env$packages          # named list of versions

# Default captures every currently loaded namespace.
names(capenv())

# the generator and whether it was seeded
set.seed(42)
capenv("stats")[c("rng_kind", "rng_seeded")]

# a package that is not installed is recorded as NA, not an error
capenv(c("stats", "notapackage123"))$packages

Description

Turns the counts a court publishes each period – cases received, cases disposed and, when it is published, cases pending – into the measures used to judge whether a court is keeping up: the clearance rate (cases disposed per case received), its rolling value over a window, the net flow into the inventory and, with the pending count, the number of periods the current inventory would take to clear at the recent pace of disposition.

Usage

caseflow(recv, disp, pend = NULL, prd = NULL, win = 12L, conflvl = 0.95)

Arguments

Details

A clearance rate above 1 means the court disposed of more cases than it received, so its inventory fell. The interval treats received and disposed as Poisson counts over the same period: given their sum, the disposed count is binomial, and the exact (Clopper-Pearson) limits for that proportion give exact limits for the ratio. The interval measures the period-to-period noise of counts this size; it does not make received and disposed the same cases, which they are not (a case disposed this month was usually received months before).

pertoclear divides the pending inventory by the mean number disposed over the last win periods (fewer at the start of the series). The published pending count need not equal the previous count plus the net flow – courts reactivate cases on bench warrants and transfer them – so penchg is reported next to net rather than derived from it.

The result carries attribute "caseflow": win, conflvl and the totals received and disposed. Errors: recv empty; a count vector that is not numeric, has a missing, negative or infinite value, or a length other than recv's; prd of another length, with missing or repeated values; win not a whole number of at least 1; conflvl outside (0, 1).

Value

A data frame of class rmbl_case_flow, one row per period, with period, received, disposed, clearance with lower and upper, net (received less disposed: positive grows the inventory), rolclr with rollo and rolup (NA until a full window is available), and either pending, penchg and pertoclear, or, without pending, cumnet (the inventory change since the first period).

References

Ostrom, B. J. and Hanson, R. A. (1999). Efficiency, Timeliness, and Quality: A New Perspective from Nine State Criminal Trial Courts. National Center for State Courts, Williamsburg, VA. (Clearance rate and time-to-disposition measures.)

Ulm, K. (1990). A simple method to calculate the confidence interval of a standardized mortality ratio. American Journal of Epidemiology, 131(2), 373-375. (Exact limits for a ratio of Poisson counts through the beta distribution.)

See Also

lifetbl, dispmix, stkflow for custody populations.

Examples

f <- read.csv(system.file(
  "extdata",
  "ocj_flow_ontario.csv",
  package = "rmoriebricklayer"
))
cf <- caseflow(f$received, f$disposed, pend = f$pending, prd = f$month)
cf
plot(cf)

# without the pending count, a quarterly rolling window and 90% intervals
cq <- caseflow(
  f$received,
  f$disposed,
  prd = f$month,
  win = 3,
  conflvl = 0.90
)
head(as.data.frame(cq)[, c("period", "clearance", "rolclr", "cumnet")])

# one value per period in every count
try(caseflow(c(100, 120, 90), c(95, 130)))

Description

Computes the BLAKE2b hash of text or bytes as lowercase hex, with any digest length from 1 to 64 bytes. With a key it becomes a message authentication code: only someone holding the key can produce or check the digest. Use it for fingerprints and keyed checks inside your own workflow; use csha256() or chmac() where other tools must reproduce the value.

Usage

cblake2b(x, key = NULL, length = 32L)

Arguments

Details

BLAKE2b (RFC 7693): a modern cryptographic hash, faster than SHA-256 in software, that takes a key natively and produces any digest length from 1 to 64 bytes.

The native key is the interesting part. A keyed BLAKE2b IS a message authentication code, with no HMAC wrapper and so no doubled hashing – cblake2b(msg, key = k) does the job of chmac() at lower cost. Use SHA-256 or HMAC where interoperability with other tools matters; use this where it does not and speed does.

Character input is hashed element by element, as the bytes of each string; NA elements give NA. A raw vector is one message. The digest length is a parameter of BLAKE2b, so a 32-byte digest is not the first half of the 64-byte digest of the same input. Errors: length outside 1 to 64; a key that is not a single string or raw vector, or is longer than 64 bytes.

Value

A character vector of lowercase hex digests – one per element for character input, length-1 for raw input.

References

Saarinen MJ, Aumasson JP (2015). The BLAKE2 Cryptographic Hash and Message Authentication Code (MAC). RFC 7693. DOI 10.17487/RFC7693

See Also

csha256(), csha512(), chmac()

Examples

# The RFC 7693 test vector for BLAKE2b-512 of "abc".
cblake2b("abc", length = 64)

# 32 bytes by default, the same width as SHA-256.
cblake2b("abc")
nchar(cblake2b("abc"))

# Any digest length, which SHA-2 cannot do.
cblake2b("abc", length = 8)

# Keyed, so it authenticates without an HMAC construction.
cblake2b("manifest", key = "secret")
cblake2b("manifest", key = "secret") ==
  cblake2b("manifest", key = "other")

# Vectorised, and raw input hashes the same bytes.
cblake2b(c("a", "b"))
identical(cblake2b("abc"), cblake2b(charToRaw("abc")))

try(cblake2b("abc", length = 65))

Description

B resamples of x, drawn with replacement and each the same length as x, with the mean of every resample returned. The spread of these replicate means estimates the standard error of the sample mean, and their quantiles give a percentile confidence interval, without assuming the data are normal (Efron and Tibshirani, 1993). The resampling uses the core's own 64-bit Mersenne Twister seeded by seed, NOT R's RNG, so a given seed reproduces the same replicates in every binding of the core and R's own random stream is left untouched.

Usage

cbootmean(x, B = 1000L, seed = 42L)

Arguments

Details

Replicate $b$ draws $n$ indices uniformly from $1, \dots, n$ with replacement (C++ std::mt19937_64 seeded with seed, then std::uniform_int_distribution, whose algorithm is the C++ library's own, so the replicates for a seed can differ between libstdc++ and libc++) and returns the mean of the drawn values, accumulated in extended precision. The replicates are a single stream: the first B replicates for a seed are the same whatever larger B is asked for later.

Edge cases: an empty x gives B values of NaN; an NA in x makes every replicate that draws it NA, so remove missing values first. B must be a whole number of at least 1 ("B must be a single whole number >= 1"); a character, factor or list x is refused ("must be numeric"). seed must be a whole number from 0 to below $2^{53}$: NA, a negative or a fractional seed is refused.

Performance: the work is $B \times n$ draws (about 0.4 s for $B = 2000$, $n = 20000$ on a laptop). The kernel polls for a Ctrl-C after about every million draws, so a long run can be interrupted; while it runs, an interactive session (or one with options(morie.progress = TRUE)) shows a spinner on stderr from 0.75 s on, and options(morie.quiet = TRUE) silences it.

Value

A numeric vector of length B: the replicate means, in the units of x.

References

Efron B, Tibshirani RJ (1993). An Introduction to the Bootstrap. Chapman and Hall.

See Also

exprspplr() for a bootstrap of a model estimate.

Examples

set.seed(1)
x <- stats::rnorm(50, mean = 5)

reps <- cbootmean(x, B = 500, seed = 7)
length(reps)

# The replicates centre on the sample mean, and their spread estimates
# the standard error.
c(sample = mean(x), bootstrap = mean(reps))
c(bootstrapse = stats::sd(reps),
  formulase = stats::sd(x) / sqrt(length(x)))

# A percentile confidence interval for the mean.
stats::quantile(reps, c(0.025, 0.975))

# Reproducible: the same seed gives the same replicates, and R's own
# random stream is not consumed.
identical(cbootmean(x, 100, seed = 1),
          cbootmean(x, 100, seed = 1))

# A longer run continues the same stream.
identical(cbootmean(x, 10, seed = 1),
          cbootmean(x, 100, seed = 1)[1:10])

# At least one replicate is needed.
try(cbootmean(x, B = 0))

Description

The covariance of every pair of columns of a numeric matrix or data frame, computed by the compiled core. Covariances describe how variables move together; the matrix is the input to principal components, Mahalanobis distances and multivariate drift checks. Column names are carried through to both dimensions of the result.

Usage

ccovm(x)

Arguments

Details

Entry $(j, k)$ is $\frac{1}{n - 1} \sum_i (x_{ij} - \bar x_j)(x_{ik} - \bar x_k)$, the n - 1 denominator of stats::cov(), with which it agrees to rounding. The diagonal holds the column variances.

Missing values are not removed: an NA in a column makes every covariance involving that column NA (pairwise or complete-case deletion is not done; use stats::cov() with ⁠use =⁠ for those). A single row gives NaN throughout. Refusals: a vector that is not a matrix ("must be a matrix or data frame") and a data frame with a non-numeric column ("every column of x must be numeric"). Logical and integer matrices are converted to double. The cost is $O(n p^2)$ for $n$ rows and $p$ columns.

Value

A symmetric ncol(x) by ncol(x) numeric matrix, with the column names of x (if any) as both row and column names; entry units are the product of the two columns' units.

See Also

ccorr() for a single Pearson correlation, cmom() for one column's variance.

Examples

X <- cbind(a = c(1, 2, 3, 4), b = c(2, 4, 7, 8), c = c(5, 3, 2, 1))
ccovm(X)

# Agrees with stats::cov().
all.equal(ccovm(X), stats::cov(X))

# The diagonal is the column variances.
all.equal(diag(ccovm(X)), apply(X, 2, stats::var))

# Data frames are accepted.
ccovm(data.frame(u = 1:5, v = c(2, 1, 4, 3, 6)))

# A plain vector is not a matrix of variables.
try(ccovm(1:5))

Description

Computes the CRC-32 checksum of character or raw input: the CRC of ITU V.42, zip and PNG (reflected polynomial 0xEDB88320, initial value and final XOR 0xFFFFFFFF). Use it as a fast first check that a large file or capsule member arrived intact, or to compare with the CRC a zip archive records for its members.

Usage

ccrc32(x)

Arguments

Details

A CRC is NOT a cryptographic digest: it detects accidental corruption – a truncated download, a flipped bit on disk – but anyone can construct a different input with the same value, so it must never be used where csha256() is meant. It is here because it is far cheaper than SHA-256 over gigabyte-scale capsule members, which makes it the right first pass when the question is only “did this file arrive intact”.

Character input is checked element by element, as its bytes, with NA giving NA; other atomic input is converted with as.character(). A raw vector gives one checksum. For files use crc32file().

Value

A numeric vector of unsigned 32-bit checksums (one per element for character input, NA for NA; length 1 for raw input). Returned as double rather than integer because values above 2^31 - 1 are not representable as an R integer.

See Also

csha256() for integrity against tampering rather than accident.

Examples

# The standard check value: CRC-32("123456789") is 0xCBF43926.
ccrc32("123456789")
ccrc32("123456789") == 0xCBF43926

# The empty input has checksum zero.
ccrc32("")

# Vectorised, and sensitive to a single changed byte.
ccrc32(c("brick", "brack"))

# Raw bytes work the same way.
ccrc32(charToRaw("123456789"))

Description

Builds the path from a leaf certificate to a trust anchor you supply, verifying each signature, each validity window, and the constraints that stop a certificate being used for something it was not issued for.

Usage

certchnvrf(
  leaf,
  trust,
  inter = list(),
  attm = Sys.time(),
  purp = NULL,
  crls = list(),
  pols = NULL,
  rvk = c("supplied", "fetch", "none")
)

Arguments

Details

attm is the point the validity windows are checked against, and it matters which one you pass. For a timestamp token the right answer is the time the token asserts: a token signed in 2020 by a certificate that expired in 2021 was validly signed, and checking it against today would reject it for no good reason. For a signature you are being asked to rely on now, pass now.

The path is built by following issuer names from the leaf through the intermediates until a certificate in trust is reached, at most ten steps (a cycle among supplied certificates cannot loop). Then every link is checked: the validity window at attm, the signature under the issuer's key (RSA PKCS#1 v1.5 with SHA-1, SHA-256, SHA-384 or SHA-512; ECDSA with SHA-256, SHA-384 or SHA-512 on P-256, P-384 and P-521), that every issuer is a CA with a path length that allows the chain, name constraints, the policy rules of RFC 5280 section 6.1, the requested purp (extended key usage of the leaf), that the anchor is self-signed, and revocation as rvk says. Each check is one row of checks; ok is TRUE only when every row passed. A failed check is a row, not an error.

Errors: trust empty; a certificate that certprs() cannot read; an unknown rvk. With revocation = "fetch" the network is used (distribution points and OCSP responders named in the certificates); otherwise nothing is fetched.

Value

A list of class bricklayer_certpath_check: ok (logical, the verdict), path (the list of certificates from the leaf to the anchor, as built), checks (a data frame with check, logical ok and detail, one row per check), attm (the time used) and policies (the outcome of policy processing). print() shows the verdict and the checks.

See Also

certprs(), tsvrf().

Examples

# The test chain shipped with the package: a self-signed P-256 root and a
# leaf it issued (a throwaway test authority; trust it for nothing). With
# your own certificates, pass their files (or certprs() results).
x509 <- function(f) {
  system.file("extdata", "x509", f, package = "rmoriebricklayer")
}
root <- certprs(x509("test-root.der"))
leaf <- certprs(x509("test-leaf.der"))

# the leaf chains to the root, checked a day after it was issued
res <- certchnvrf(leaf, trust = root, attm = leaf$not_before + 86400)
res$ok
res$checks[, c("check", "ok")]

# after the leaf expired, the same chain fails on validity
certchnvrf(leaf, trust = root, attm = leaf$not_after + 86400)$ok

# a purpose the leaf was not issued for
r <- certchnvrf(
  leaf,
  trust = root,
  attm = leaf$not_before + 86400,
  purp = "timeStamping"
)
r$checks[!r$checks$ok, ]

# a chain must end at something you trust
try(certchnvrf(leaf, trust = list()))


# with revocation looked up on the network
certchnvrf("tsa.crt", trust = "ca.crt", purp = "timeStamping",
                  rvk = "fetch")

Description

Reads the fields a verifier needs: who issued it, who it is for, when it is valid, what key it carries, and what that key is allowed to do. Use it to inspect the certificate of a timestamping authority or a signer before you decide to trust it, or as the first step of certchnvrf(). The DER decoder is the package's own; no system library or network is used.

Usage

certprs(cert)

Arguments

Details