--- title: "Performing a typed analytical plan" output: rmarkdown::html_vignette vignette: > %\VignetteIndexEntry{Performing a typed analytical plan} %\VignetteEngine{knitr::rmarkdown} %\VignetteEncoding{UTF-8} --- ```{r setup, include = FALSE} knitr::opts_chunk$set(collapse = TRUE, comment = "#>") ``` ```{r library} library(conductoR) ``` ## The question Does nitrogen rate affect maize yield? `agridat::hernandez.nitrogen` holds 136 plots from five Minnesota sites, four replicates each, at nitrogen rates from 0 to 268 kg/ha. Two methods answer the question independently: the slope of a linear regression of yield on rate, and Kendall's rank correlation. A plan puts a power-aware gate on the rank lens, reconciles the two lenses, and turns the consensus into a recommendation. ```{r data} data("hernandez.nitrogen", package = "agridat") d <- hernandez.nitrogen table(d$site, d$rep) ``` ## The tools A tool is a function of `(inputs, data, node)` that returns an `orchestraManifest::orchestra_manifest`. `lens_to_manifest()` lifts a yes-or-no verdict into one; `decision_to_manifest()` does the same for a recommendation. The tools are registered by name so the plan can be written to a file. ```{r tools} ols_lens <- function(inputs, data, node) { fit <- lm(yield ~ nitro, data = data) co <- summary(fit)$coefficients["nitro", ] p <- unname(co["Pr(>|t|)"]) lens_to_manifest( headline = if (p < 0.05) "effect" else "no_effect", effect_present = p < 0.05, emitter = "ols", method = "ols:slope_t", metrics = list(p_value = p, slope = unname(co["Estimate"]), effect_mean = unname(co["Estimate"]) * diff(range(data$nitro)), noise_sd = sigma(fit), n = nrow(data))) } rank_lens <- function(inputs, data, node) { ct <- suppressWarnings(cor.test(data$nitro, data$yield, method = "kendall")) lens_to_manifest( headline = if (ct$p.value < 0.05) "effect" else "no_effect", effect_present = ct$p.value < 0.05, emitter = "rank", method = "rank:kendall", metrics = list(p_value = ct$p.value, tau = unname(ct$estimate), n = nrow(data))) } recommend <- function(inputs, data, node) { consensus <- inputs$consensus@summary$headline action <- switch(consensus, causal_effect_corroborated = "apply_nitrogen", no_effect_corroborated = "withhold_nitrogen", "insufficient_evidence") decision_to_manifest(action, emitter = "agronomy", method = "agronomy:rule", metrics = list(consensus = consensus), consumed = list(inputs$consensus)) } register_tool("ols_lens", ols_lens) register_tool("rank_lens", rank_lens) register_tool("recommend", recommend) ``` ## The plan Five nodes. Each edge names the inferential target it expects from its source; the gate reads the rank lens and abstains when a null verdict comes from fewer than 25 plots or from a design that could not have resolved the implied effect. ```{r plan} score <- oap_plan( id = "nitrogen_response", meta = list(question = "Does nitrogen rate affect maize yield?"), nodes = list( node("lens_ols", "ols_lens", emits = "treatment_effects", role = "regression"), node("lens_rank", "rank_lens", emits = "treatment_effects", role = "rank"), node("gate", operator = "abstain_gate", emits = "decisions", consumes = list(edge("lens_rank", "treatment_effects")), params = list(n_min = 25L, z_floor = 2, mode = "flag")), node("consensus", operator = "triangulate", emits = "decisions", consumes = list(edge("lens_ols", "treatment_effects"), edge("lens_rank", "treatment_effects"), edge("gate", "decisions"))), node("decision", "recommend", emits = "decisions", role = "agronomy", consumes = list(edge("consensus", "decisions"))) )) plan_order(score) ``` ## Scenario A: the full trial ```{r scenario_a} run_a <- conduct(score, data = d) run_a$terminal@summary$headline run_a$bus$consensus@summary$metrics[c("agree", "n_lenses")] ``` Both lenses find the effect, the gate passes, and the consensus is corroborated. ## Scenario B: one site, one replicate, three rates ```{r scenario_b} small <- d[d$site == "S1" & d$rep == "R1" & d$nitro <= 67.2, ] nrow(small) run_b <- conduct(score, data = small) run_b$bus$gate@summary$abstain_reason run_b$terminal@summary$headline ``` With three plots the rank test cannot reject, and the gate refuses to let that null through. In `flag` mode the consensus abstains and the recommendation is `insufficient_evidence`. Rerouting instead drops the gated lens and reconciles what survives. ```{r scenario_b_reroute} gate <- plan_node(score, "gate") gate@params$mode <- "reroute" rerouted <- score rerouted@nodes[[3]] <- gate run_b2 <- conduct(rerouted, data = small, verbose = FALSE) run_b2$bus$consensus@summary$metrics[c("rerouted", "n_lenses")] run_b2$terminal@summary$headline ``` One lens on three plots is not a corroborated effect, so the recommendation is still to wait for evidence. ## The plan as a file Because every compute node names its tool, the plan writes to YAML and reads back to the same run. ```{r yaml} f <- tempfile(fileext = ".oap.yaml") plan_to_yaml(score, f) cat(readLines(f, n = 12), sep = "\n") again <- plan_from_yaml(f) run_c <- conduct(again, data = d, verbose = FALSE) identical(run_c$terminal@run_id, run_a$terminal@run_id) ``` ## Writing a plan A plan can also be written directly as a file, by a person or by an AI model. Version 0.2 of the plan language adds four optional fields that let a plan act on its own results: | Field | What it does | |---|---| | `when` | runs the node only if a condition on upstream results is true | | `foreach` | runs the node once per item, from a list or from a column of an upstream result | | `on_refusal` | says what an abstaining result does: `continue`, `halt-branch` or `stop` | | `operator: select` | keeps the rows of an upstream table that meet a rule | The plan below asks which interactions in a two-way model of yield hold up. The first tool fits `yield ~ nitro * site` and returns one row per model term. A `select` node keeps the interaction terms with p < 0.05, and a check runs once for each kept term, but only if there is at least one. ```{r term_tools} term_table <- function(inputs, data, node) { a <- anova(lm(yield ~ nitro * site, data = data)) terms <- rownames(a)[rownames(a) != "Residuals"] tab <- data.frame(term = terms, order = lengths(strsplit(terms, ":", fixed = TRUE)), p_value = a[terms, "Pr(>F)"]) summ <- orchestraManifest::manifest_summary("terms", metrics = list(n = nrow(data))) orchestraManifest::orchestra_manifest( emitter_package = "example", emitter_version = "0.1.0", inferential_target = "parameters", run_id = "terms-1", method = "lm:anova", params = tab, summary = summ, timestamp = Sys.time(), data_hash = orchestraManifest::manifest_data_hash(tab, NULL, NULL, NULL, NA_integer_, summ)) } ``` The check is a stand-in for an external validation: it refits the model on replicates R1 and R2, then on R3 and R4, and asks whether the term is significant in both halves. When a half has too few plots it abstains with `insufficient_evidence`, which `on_refusal: continue` records and passes on. ```{r check_tool} check_term <- function(inputs, data, node) { term <- node@params$item p <- vapply(list(c("R1", "R2"), c("R3", "R4")), function(half) { a <- anova(lm(yield ~ nitro * site, data = data[data$rep %in% half, ])) a[term, "Pr(>F)"] }, numeric(1)) verdict <- if (anyNA(p)) "insufficient_evidence" else if (all(p < 0.05)) "replicated" else "not_replicated" decision_to_manifest(verdict, emitter = "split_half", method = paste0("split_half:", term), metrics = list(p_first = p[1], p_second = p[2])) } register_tool("term_table", term_table) register_tool("check_term", check_term) ``` The plan as YAML. Every string is in quotes: YAML 1.1 reads an unquoted `yes`, `no`, `on` or `off` as a logical, so conductoR refuses those words unquoted rather than guess. ```{r plan_v02} plan_text <- ' id: "interactions" plan_version: "0.2" meta: question: "Which nitrogen-by-site interactions hold up in both halves of the trial?" nodes: - id: "terms" tool: "term_table" emits: "parameters" - id: "flagged" operator: "select" emits: "parameters" consumes: - from: "terms" contract: "parameters" params: rule: "order > 1 & p_value < 0.05" - id: "check" tool: "check_term" emits: "decisions" consumes: - from: "flagged" contract: "parameters" when: "flagged.n_rows > 0" foreach: from: "flagged" column: "term" foreach_max: 10 on_refusal: "continue" ' interactions <- plan_from_yaml(plan_text) run_i <- conduct(interactions, data = d) run_i$conditions run_i$items[, c("item", "status")] run_i$bus$check@outputs[["check[nitro:site]"]]@summary$headline ``` One interaction, `nitro:site`, meets the rule on all 136 plots, and the check finds it significant in the first half but not in the second. The run record lists the condition, its value and every item with its status, so a blinded run and the unblinded run that follows it can be compared item by item. Conditions are read by conductoR, never run as R code. They may refer to `.headline`, `.abstained`, `.n_rows` and `.metrics.`, and use plain numbers, quoted strings, comparisons, `%in%` with `c(...)`, `&`, `|`, `!` and brackets. Anything else is refused when the plan is read: ```{r refused} res <- try(plan_from_yaml(sub("flagged.n_rows > 0", "system(cmd) == 0", plan_text, fixed = TRUE)), silent = TRUE) cat(conditionMessage(attr(res, "condition"))) ``` ### The schema and the repair loop `plan_schema()` returns a JSON Schema (draft 2020-12) of the language. An editor can check a plan against it as it is typed, and an AI model can use it to write a plan as structured output. The loop is: write the plan as JSON against the schema, check it with `plan_validate()`, and fix every problem the table lists. `plan_validate()` reports all the problems at once, each with the field, the node, the rule broken and the fix. ```{r schema} schema <- jsonlite::fromJSON(plan_schema(), simplifyVector = FALSE) names(schema$`$defs`$node$properties) draft <- '{ "id": "interactions", "plan_version": 0.2, "nodes": [ {"id": "terms", "tool": "term_table", "emits": "parameters"}, {"id": "flagged", "operator": "select", "emits": "parameters", "consumes": [{"from": "terms", "contract": "parameters"}], "params": {"rule": "order > 1 & p_value < 5e-2"}}, {"id": "check", "tool": "check_term", "consumes": [{"from": "flagged", "contract": "parameters"}], "foreach": {"from": "flagged", "column": "term"}, "on_refuse": "continue"} ] }' plan_validate(draft)[, c("node", "field", "fix")] ``` Each row says what to change. With the version quoted, the number written plainly, `emits` declared and the field name corrected, the plan is valid, and `plan_to_yaml()` writes the YAML a person reads. ```{r schema_fixed} fixed <- draft fixed <- sub('"plan_version": 0.2', '"plan_version": "0.2"', fixed, fixed = TRUE) fixed <- sub("5e-2", "0.05", fixed, fixed = TRUE) fixed <- sub('"tool": "check_term",', '"tool": "check_term", "emits": "decisions",', fixed, fixed = TRUE) fixed <- sub('"on_refuse"', '"on_refusal"', fixed, fixed = TRUE) nrow(plan_validate(fixed)) cat(plan_to_yaml(plan_from_json(fixed))) ``` ## A mistyped edge An edge that declares the wrong contract is refused before the downstream node runs. Here the decision node is told to expect `parameters` from the consensus, which emits `decisions`. ```{r violation} bad <- oap_plan(id = "mistyped", nodes = list( node("lens_ols", "ols_lens", emits = "treatment_effects"), node("decision", "recommend", emits = "decisions", consumes = list(edge("lens_ols", "parameters"))) )) res <- try(conduct(bad, data = d, verbose = FALSE), silent = TRUE) cat(conditionMessage(attr(res, "condition"))) ``` ## Provenance and a targets pipeline Giving `conduct()` an `out_dir` writes `run.rds`, `run_card.md` and `provenance.md`, the last carrying the chain of payload hashes in execution order. `run_provenance()` returns the same table in the session. ```{r provenance} run_provenance(run_a)[, c("node", "operator", "emits", "abstained")] ``` The same plan compiles to a `_targets.R` in which each node is one target calling `conduct_step()`, so a long-running plan inherits caching from the 'targets' package. ```{r targets, eval = FALSE} compile_to_targets(score, "_targets.R", preamble = c( "library(conductoR)", "source('tools.R')", ".score <- plan_from_yaml('nitrogen_response.oap.yaml')", ".data <- agridat::hernandez.nitrogen")) targets::tar_make() targets::tar_read_raw(terminal_target(score)) ```