What difference can a sponsored place on a bioinformatics course make? For some researchers, it means access to training they would otherwise struggle to afford. For others, it is an opportunity to explore a new field, work through unfamiliar analytical methods, or find new ways to approach an existing research question.

Through the Bioconductor–Physalia collaboration, Physalia offers two free places on eligible courses to researchers based in low- and middle-income countries (LMICs). The Bioconductor Community Manager coordinates the application process and works with the Training Committee co-chairs on selecting participants, while Bioconductor helps promote the courses through its community channels. Between November 2025 and October 2026, the programme provided 18 sponsored places across nine courses, covering topics ranging from RNA-seq and proteomics to microbiome analysis, epigenomics and network biology.

The aim is to make specialist training more accessible to researchers who might otherwise miss out on these opportunities. We asked participants what they found useful about their courses and how they planned to use what they had learned. Their responses offer a glimpse into the different ways the training connects with their work.

© 2026 Bioconductor. Content is published under Creative Commons CC-BY-4.0 License for the text and BSD 3-Clause License for any code. | R-Bloggers

I’ve been giving permute a bit of attention recently. Two additions should make it easier to set up permutation tests: an option for ordered set partitions, and formulas for specifying plots and blocks. That’s a rather grand name for the first addition, so let’s start with what it does.

These changes are in the development version. To follow along, install it from my R-universe:

install.packages("permute", repos = c(
  "https://gavinsimpson.r-universe.dev",
  "https://cloud.r-project.org"
))

library("permute")
data(jackal, package = "permute")

table(jackal$Sex)

  Male Female 
    10     10 

numPerms(nrow(jackal), control = how(within = Within(type = "free")))

[1] 2.432902e+18

ctrl <- how(
  plots = Plots(strata = ~ Sex, type = "partition"),
  data = jackal
)

numPerms(nrow(jackal), control = ctrl)

[1] 184756

set.seed(42)
p <- shuffle(nrow(jackal), control = ctrl)

table(jackal$Sex[p])

  Male Female 
    10     10 

fit <- lm(Length ~ Sex, data = jackal)
coef(fit)

(Intercept)   SexFemale 
      113.4        -4.8 

set.seed(42)
perms <- shuffleSet(nrow(jackal), nset = 999, control = ctrl)

observed <- unname(coef(fit)[2])
null <- apply(perms, 1, function(i) {
  permuted <- jackal
  permuted$Sex <- jackal$Sex[i]
  unname(coef(lm(Length ~ Sex, data = permuted))[2])
})

hist(null, breaks = "FD", main = "",
     xlab = "Female minus male mean mandible length (mm)",
     col = "grey85", border = "white")
abline(v = observed, col = "firebrick", lwd = 2)

tol <- sqrt(.Machine$double.eps)
p_value <- (1 + sum(abs(null) >= abs(observed) - tol)) /
  (nrow(perms) + 1)

c(difference_mm = observed, p_value = p_value)

difference_mm       p_value 
       -4.800         0.005 

library("vegan")
data(dune, package = "vegan")
data(dune.env, package = "vegan")

formula_ctrl <- how(
  plots = Plots(strata = ~ Management / Use),
  blocks = ~ Management,
  data = dune.env
)

set.seed(42)
i <- shuffle(nrow(dune.env), control = formula_ctrl)

all(getBlocks(formula_ctrl)[i] == getBlocks(formula_ctrl))

[1] TRUE

all(getStrata(formula_ctrl)[i] == getStrata(formula_ctrl))

[1] TRUE

env <- dune.env
env$group <- getStrata(formula_ctrl)
ord <- rda(dune ~ A1 + Condition(group), data = env)

set.seed(42)
anova(ord, permutations = formula_ctrl)

Permutation test for rda under reduced model
Blocks:  Management 
Plots: Management/Use, plot permutation: none
Permutation: free
Number of permutations: 199

Model: rda(formula = dune ~ A1 + Condition(group), data = env)
         Df Variance      F Pr(>F)
Model     1   4.3492 1.3265  0.225
Residual  8  26.2298              

plot_groups <- Plots(strata = ~ Management / Use, data = dune.env)

Claude Code (and many other agents) has been getting so good at ecological analysis that I find I’m increasingly tempted not to directly review the code it writes. In fact, it finds more mistakes in my code than I find in its code these days.

So I need a way to check the validity of the code it writes. Here’s the strategy for this I’ve been developing. Its a work in progress.

Rico Mariani, a long-time performance engineer at Microsoft, coined the idea of the pit of success. A system has a pit of success when the natural, lazy way of using it leads to good outcomes. You don’t need heroics or perfect memory. You fall into the right result, and climbing out to do something wrong takes deliberate effort.

I came across this framing in a talk that I was recently recommended and I further recommend it to anyone interested in the topic: Functional architecture – The pits of success – Mark Seemann. It put a name to something I had been circling around for years.

I had been arguing a version of the same point long before I knew the term. In my 2022 post titled Functional programming explains why containerization is needed for reproducibility, I argued that a pipeline is only as reproducible as the functions it is built from, and that even a pure function hides inputs such as the R version, the package versions, and the operating system. Containers are one way to make those hidden inputs explicit. In another 2022 post on learnings from functional programming, I made a related case for writing functions that fail early and never reach into the global environment. Both posts were, in a sense, attempts to build a pit of success out of habits and tooling. Neither was framed that way, and I think the framing is more useful.

“This original, thought-provoking, and transforming, book is much much more than an implementation manual for Bayesian Data Analysis, even though it shares almost the same perspective. (The first sentence of the book states that the authors’ `conceptions of statistical practice, and of Bayesian statistics, have changed over the years’.) By providing a modus vivendi for undertaking Bayesian modelling from scratch in realistic settings where models are not magicked out of the blue, the authors explicit and rationalise the many steps required by such a bottom-up modelling protocol (`not a checklist, not a cookbook’, and not a flowchart!) in real situations. The contents read very well and very smoothly, with a seamless conjunction of intuition, modelling advices, computational details, and comparison tools. While unsurprisingly Bayesian, the perspective adopted therein remains both open and inclusive, with a welcome humility about the limitations and challenges of Bayesian workflows. This book should thus appeal to and profit a wide variety of readers, as providing guidance through an extensive collection of highly detailed examples, with shared code and exercises.”

This book proposes a modus vivendi for Bayesian modelling in applied, realistic Bayesian analysis, where models are not magicked out of the blue. It thus emphases iterative model building, model checking, computational troubleshooting, and simulated-data experimentation, filling a gap that looks glaring in retrospect. It particularly targets users and developers of Stan, with code excerpts in R and Stan. It consists of four parts:

1. background on Bayesian methods and computational tools;
2. the Bayesian workflow proper, namely building a statistical model from its components, together with its assessment tools;
3. the computational aspects of fitting models, diagnosing convergence and assessing calibration;
4. case studies.

I was eagerly waiting for the book, as I knew Andrew, Aki, and Richard had been working on it for a few years. (The quote above is the blurb I wrote upon request from the publisher.)

The tenets of BaWoFlo—if I may resort to this acronym!—are (i) fitting multiple models, (ii) applying methods repeatedly, and (iii) resorting to simulated-data experiments, which should not come as a surprise to readers of BDA. As noted in the introduction, the protocol exposed therein can also benefit non-Bayesian experimenters. This agrees with the highly moderate, “M-open”, agnostic approach to Bayesianism adopted by the authors (“there is no safe haven”). I also welcome and share their humble perspective about the limitations and challenges of Bayesian workflows.

Examples are treated in full detail, with successive modelling and computational choices profusely commented, which is a big plus for such a practical book. This starts as early as Chapter 4, with a multiple-choice exam example. Indeed, there cannot be general principles or a generic theory that would make the approach foolproof. See, e.g., “A data model is not just a ‘likelihood’” (p.70), as when the data model is not fully generative. I very much liked the section on choosing priors (5.6), and the very rich graphs (see, e.g., Chapter 8) for assessing the impact of prior and likelihood, as well as for predictive checks. In coherent continuation of the authors’ earlier work, the book advocates LOO methods and model stacking rather than model averaging. (With a surprisingly anti-Ockham perspective in Section 9.7.)

The MCMC coverage is unsurprising, with \(\hat R\) at the forefront. Chapter 12, on using fast experiments to detect fitting or computational issues, is very nice. The book builds on the immense corpus of work achieved by the authors over the decades (for the most senior ones!). By contrast, the chapter on approximate solutions (13) is way too short, and the same goes for those on calibration and software development.

The book is very US-centric, unsurprisingly given Andrew’s focus on political science. Some sections are reminiscent of Andrew’s blog entries (or the opposite). The (football) World Cup example was initiated when Andrew was in France, during the 2014 World Cup, and as a result (?) the names of the teams are in French! One chapter also reanalyses the birthdate data displayed on the cover of BDA.

Mileage varies on the applied chapters, depending on the example. A dog chapter is followed by a cat chapter! Not that the (stat)dog experiment was in any way enjoyable, especially for the dogs. Maybe the cats were running it! And then come chapters on roaches and sharks. There is also a frightening flowchart (Fig. 2.1)! And the book ends with an appendix on going through BDA to better understand BaWoFlo

[The usual disclaimer applies, namely that this review is likely to appear later in CHANCE, in my book reviews column.] 

On 15-17 September 2026, members of the Bioconductor, Galaxy and nf-core communities came together at The Open University in Milton Keynes for the BioFAIR Data to Discovery and BioFAIR Collaboration Fest (CoFest) joint event.

The three-day event combined hands-on training with collaborative work on research software interoperability. It was jointly organised by BioFAIR Fellow Marisa Loach and BioFAIR Pathfinder project lead and Institute for Research Software 2026 Fellow Kevin Rue-Albrecht.

The CoFest brought together contributors based in the UK, Ireland, Germany, India and the USA, with both in-person and remote participation.

Marisa led the Data to Discovery training, introducing participants to single-cell analysis using the European instance of the Galaxy platform (https://usegalaxy.eu/), while Kevin led the hybrid CoFest focused on making methods developed as Bioconductor R packages more readily accessible through infrastructures such as Galaxy and nf-core.

A Galaxy community recap of the event also covers the training and CoFest from the Galaxy perspective.

© 2026 Bioconductor. Content is published under Creative Commons CC-BY-4.0 License for the text and BSD 3-Clause License for any code. | R-Bloggers

Before the main content: * With all the funding cuts to PhD studies, I’d appreciate if you can donate here.

I am using Omarchy, a Linux distribution based on Arch for those that want efficiency, and it is similar to Manjaro minus the desktop effects that use considerable resources.

At Surrey and probably other universities in the UK, we require to use Microsoft Authenticator to log to the VPN.

Today I had to connect to Eureka2, a powerful computing cluster, using my personal laptop. The official guide Working off campus with GlobalProtect instructs to download Global Protect from the Surrey App Store, a service that has been retired.

I installed Open Connect which is a FOSS alternative:

sudo yay -S openconnect

I had to adjust my email settings here, and verify that the default sign-in method is ‘Microsoft Authenticator – notification’.

Then I was able to connect using:

sudo openconnect --protocol=gp --usergroup=gateway vpn.surrey.ac.uk

After a few tries, what works it to use the full username (not email name) with ‘@surrey.ac.uk’. What I mean is that my university account is something like ab12345 when my email firstname.lastname, and so I used ab12345@surrey.ac.uk. For the password I had to use the password that I use for Surrey Learn (not the email password) followed by the authenticator timed PIN number.

I was finally able to see:

ESP session established with server
ESP tunnel connected; exiting HTTPS mainloop.

And finally I could type this in another terminal:

ssh mv00622@eureka2.surrey.ac.uk

Which shows:

    +++++++++++++++++++++++++++++++++++++++++
      Welcome to Surrey's Eureka2 HPC Cluster
             Based on Rocky 8
    +++++++++++++++++++++++++++++++++++++++++

Documentation for HPC users is at below site:
https://docs.pages.surrey.ac.uk/research_computing/

I hope this helps. I had to spend over half an hour figuring out the steps by reading different sources and I hope this helps other people.

There has been an explosion of package/environment/platform managers in the R ecosystem recently,¹ and R users are really reaping the rewards of so many ripe and riveting tools. But, I’ve found it rough² to keep track of them all and understand what they are and how they are different!

I suspect I am not the only one. My colleague Libby Heeren and I recently held a Data Science Lab with Eric Nantz to showcase {rix} and Nix (which was awesome! Keep an eye out on the Posit YouTube channel for when the recording is published), and the first three or so questions from the audience were all about how {rix} compares to X ({renv}, rv, uvr), and so I thought it’d be helpful to write a short summary of five tools that I’ve heard about: {renv}, {rix}, rv, uvr, and ir.³

While they’re all about tracking or isolating dependencies, they operate at a different “levels” from each other. TL;DR:

• {rix} snapshots the entire system environment
• {renv}, rv, and uvr manage R package dependencies within a single R project
• ir handles dependencies for a single R script rather than a project

I don’t want to/can’t pass judgement on which is the best or which I recommend, because I honestly haven’t had the opportunity to explore all of them fully (a very diplomatic way of saying that I haven’t yet used all of them in a real-world context), but I would love to hear everybody else’s opinions on which they think is the best tool because I imagine people have them!

It’s an unsatisfactory proxy for quality/experience, but here is the GitHub star trajectory of each:

(Mind the recency bias).

With that, let’s begin, going chronologically by release!

Terminal

install.packages("renv")

Console

renv::init()

Console

renv::snapshot()

renv.lock

{
  "R": {
    "Version": "4.4.2",
    "Repositories": [
      {
        "Name": "CRAN",
        "URL": "https://cloud.r-project.org"
      }
    ]
  },
  "Packages": {
    "dplyr": {
      "Package": "dplyr",
      "Version": "1.1.4",
      "Source": "Repository",
      "Repository": "CRAN",
      "Hash": "fedd9d00c2944ff00a0e2696ccf048ec"
    },
    "ggplot2": {
      "Package": "ggplot2",
      "Version": "3.5.1",
      "Source": "Repository",
      "Repository": "CRAN",
      "Hash": "44c6a2f8202d5b7e6c72a95ca7ea1795"
    }
  }
}

Console

renv::restore()

Terminal

install.packages("rix")

generate_env.R

library(rix)

rix(
  r_ver = "4.4.2",
  r_pkgs = c("dplyr", "ggplot2"),
  project_path = "."
)

Terminal

nix-build

Terminal

brew install rv-r

Terminal

rv init hello_world

hello_world/
  rv/
  library/
  scripts/
    activate.R
    rvr.R
  .gitignore
  .Rprofile
  rproject.toml

rproject.toml

r_version = "4.4.2"

# A list of packages to install and any additional configuration
# Examples:
    # "dplyr",
    # {name = "dplyr", repository = "CRAN"},
    # {name = "dplyr", git = "https://github.com/tidyverse/dplyr.git", tag = "v1.1.4"},
dependencies = [
    "dplyr",
    "ggplot2"
]

Terminal

rv sync

rv.lock

[[packages]]
name = "dplyr"
version = "1.1.4"
repository = "CRAN"
source = "Repository"
hash = "fedd9d00c2944ff00a0e2696ccf048ec"

[[packages]]
name = "ggplot2"
version = "3.5.1"
repository = "CRAN"
source = "Repository"
hash = "44c6a2f8202d5b7e6c72a95ca7ea1795"

Terminal

curl -fsSL https://raw.githubusercontent.com/nbafrank/uvr/main/install.sh | sh

Terminal

mkdir my-project && cd my-project
uvr init --r-version ">=4.4.2"

Terminal

uvr add dplyr ggplot2

Terminal

uvr sync

uvr.lock

r-version = "4.4.2"

[[packages]]
name = "dplyr"
version = "1.1.4"
source = "cran"
hash = "fedd9d00c2944ff00a0e2696ccf048ec"

[[packages]]
name = "ggplot2"
version = "3.5.1"
source = "cran"
hash = "44c6a2f8202d5b7e6c72a95ca7ea1795"

Terminal

uvr run analysis.R

Terminal

uvr sync --frozen

Terminal

curl -fsSL https://raw.githubusercontent.com/r-lib/ir/main/scripts/install.sh | sh

script.R

#!/usr/bin/env -S ir run
#| packages:
#|   - dplyr>=1.0
#|   - ggplot2@3.5.1
#| isolated: true
#| exclude-newer: "2025-05-15"

library(dplyr)
library(ggplot2)

mtcars |> count(cyl, gear) |> ggplot(aes(x = cyl, y = gear)) + geom_col()

report.qmd

---
title: My report
ir:
  packages:
    - dplyr>=1.0
    - ggplot2@3.5.1
  isolated: true
  exclude-newer: 2025-05-15
---

Terminal

ir run script.R
ir render report.qmd

After 20 years of service the R-Forge package hosting platform will be retiring soon. Codeberg is a modern alternative that shares many of the ideas and values that R-Forge was built on. A tutorial for transitioning to Codeberg (or other Git-based platforms) is provided.

You say goodbye

R-Forge was created 20 years ago after the 2nd useR! conference in 2006. Kurt Hornik had the idea that “R should have its own SourceForge!” (note that GitHub did not yet exist at the time). However, the goal was that this code sharing platform should be created by the R community for the R community. The implementation was led by Stefan Theussl in the following months using GForge/FusionForge, a fork of the software behind SourceForge. Hosting and funding was provided by WU Wirtschaftsuniversität Wien and The R Project.

R-Forge offers a Subversion (SVN) repository for collaborative development of R packages along with various other services that were state of the art at the time: daily package builds and checks, forums, mailing lists, bug trackers, web pages. A short introduction was published in the inaugural issue of The R Journal: Theußl & Zeileis (2009). “Collaborative Software Development Using R-Forge.” The R Journal, 1(1):9-14. doi:10.32614/RJ-2009-007

As the FusionForge system is no longer maintained, R-Forge will be retired in the next months, switching to a read-only site. The system was somewhat dated anyways: The website has a certain “00s” vibe (see the screenshot of my profile page below) and it lacks services that became standard over the last two decades like issues, pull requests, and continuous integration.

I say hello

Codeberg is a non-profit, community-led effort that provides services to free and open-source projects, such as Git hosting (using Forgejo), pages, continuous integration, etc. It is run by the non-profit organization Codeberg e.V. based in Berlin, Germany, which contributors can (but do not have to) become a member of. (I joined this summer.)

Thus, Codeberg shares many of the ideas and values, that R-Forge was built on, but provides state-of-the-art open-source tools for collaborative development of R packages (and beyond). For a first impression of the look and feel, see the screenshot of my profile page below.

I say go go go

For transitioning from R-Forge to Codeberg – or some other Git-based code hosting platform – Reto Stauffer has written up a step-by-step guide. The most important steps are also listed below.

1. SVN checkout from R-Forge

The starting point for the transition is the SVN repository from R-Forge that you probably have already checked out. For illustration we are using the c403 R package that is one of the packages contained in the exams project on R-Forge.

svn checkout svn://R-Forge.R-project.org/svnroot/exams/

2. List of R-Forge contributors

To find all SVN usernames who ever contributed to the SVN on R-Forge repository (in this example, to the folder exams/pkg/c403) the command svn log -q can be used. Using the following shell command a text file with the unique names can be created. Note that exams/pkg/c403 needs to be replaced for other projects and the output file name can be adjusted too:

(cd exams/pkg/c403 && svn log -q | grep -e '^r' | awk 'BEGIN { FS = "|" } ; { print $2 }' | sort | uniq) > \
  authors_c403.txt

This creates the text file authors_c403.txt which will look something like:

  azeileis
  retos

This file needs to be modified by adding the full names and e-mail addresses of the contributors. Important: The leading spaces need to be removed in each line for the subsequent commands to work.

azeileis = Achim Zeileis <Achim.Zeileis@R-project.org>
retos = Reto Stauffer <Reto.Stauffer@uibk.ac.at>

Provided the profile pages on R-Forge are up to date for the users, it is also possible to scrape this information, see the optional step 2*. below.

3. Clone a Git version of the SVN

Using git-svn a Git version of the SVN can be easily created. For this step it is crucial that all users are listed in the authors file created in the previous step, otherwise git-svn will fail.

git svn clone svn://R-Forge.R-project.org/svnroot/exams/pkg/c403 \
    --authors-file=authors_c403.txt --no-metadata \
    git-c403

Again, adjust the URL of the SVN, the authors file name, and the name of the git repository (git-c403) as needed.

4. Optional: Filter files in the Git repository

This step is not needed for the c403 example but is particularly useful when migrating a full R-Forge repository (rather than just a pkg sub-folder). In that case it might be useful to filter out certain files or directories such as the initial README file auto-created by R-Forge or the www sub-directory (typically just containing original PHP file or a auto-generated web page from pkgdown or altdoc).

Such files can be filtered, i.e., erased completely from the repository history, using git-filter-repo:

git-filter-repo --force --path 'README' --invert-paths
git-filter-repo --force --path 'www' --invert-paths

5. Import Git repository on Codeberg

First, create a new project on Codeberg, e.g., https://codeberg.org/zeileis/c403. The project should be completely empty without any initial commit.

Then, all left to do then is to tell the Git repository where its new home is (again adjust all names correspondingly):

cd git-c403
git remote add origin ssh://git@codeberg.org/zeileis/c403.git
git push --set-upstream origin main

Analogous commands can be used for setting the origin to a project on GitHub, GitLab, etc.

All history will be immediately available in the new repository. If the authors have accounts on Codeberg, their profiles will be linked based on the e-mail address(es) associated with their Codeberg account.

6. Package builds and checks via R-universe

To automatically build and check the development versions of the R package on Codeberg, the most convenient solution is to use the R-universe service developed by Jeroen Ooms and provided by rOpenSci.

In short, to make use of it, you need a GitHub account and create a repository like zeileis.r-universe.dev. This should contain a packages.json file with the names and repository URLs of one or more R packages, e.g.,

[
    {
        "package": "c403",
        "url": "https://codeberg.org/zeileis/c403"
    }
]

As a final step, the R-universe app has to be installed on the same GitHub account. See the R-universe set up documentation for more details.

7. Codeberg pages

If you want to add a project page for your package, e.g., generated via altdoc or pkgdown, then https://codeberg.page/ can be used.

The c403 package does not have a page, yet, but for an example see https://zeileis.codeberg.page/topmodels/ which is generated from the topmodels R package using altdoc in combination with quarto.

All that needs to be done to create the page manually, is to create an empty branch pages via git switch --orphan pages, adding the static site files, and pushing them to the repository.

Alternatively, it is also possible to use Forgejo actions to do this automatically.

8. Update descriptions and links

Finally, inform users where the sources of the package can be found, i.e., update the DESCRIPTION (especially the URL and BugReports fields) and in the NEWS. Additionally, a README.md file would be expected in a project while on R-Forge this was often not used (because R-Forge would be aware of the DESCRIPTION instead). Finally, it is recommended to add a link from the old R-Forge project to the new Codeberg project (possibly also in the www folder). To assure that none of the contributors can make modifications in the old R-Forge SVN, the permissions on R-Forge can be switched to read-only (typically for the Senior Developer role).

Optional addition: 2.* Scrape author names and e-mail addresses

In step 2, after the file with the author names (authors_c403.txt in the example above) has been created, one can try to automate the addition of proper names and e-mail addresses. This is particularly useful if there are more than just a handful of users. It works when the R-Forge user pages list the proper names and the (current) e-mail addresses.

The function below can be called in R as

rforge_authors("c403")

It first reads the authors_c403.txt file, looks up all users on R-Forge (using the rvest package), and then updates the authors_c403.txt file accordingly.

rforge_authors <- function(x) {
  txt <- sprintf("authors_%s.txt", x)
  if (file.exists(txt)) {
    a <- readLines(txt)
    a <- gsub(" ", "", a)
  } else {
    stop(sprintf("Please use 'svn log' to set up an '%s' file with all user names", txt))
  }

  for (i in seq_along(a)) {  
    doc <- a[i] |>
      sprintf(fmt = "https://r-forge.r-project.org/users/%s/") |>
      rvest::read_html()
    n <- doc |>
      rvest::html_table() |>
      lapply(as.data.frame) |>
      lapply(subset, subset = X1 == "Real Name:", select = "X2", drop = TRUE) |>
      unlist() |>
      unique()
    if (length(n) < 1L) n <- "FIXME"
    e <- doc |>
      rvest::html_table() |>
      lapply(as.data.frame) |>
      lapply(subset, subset = X1 == "Email Address:", select = "X2", drop = TRUE) |>
      unlist() |>
      unique() |>
      gsub(pattern = " @nospam@ ", replacement = "@", fixed = TRUE)
    if (length(e) < 1L) e <- "FIXME"
    a[i] <- sprintf("%s = %s <%s>", a[i], n, e)
  }
  writeLines(a, txt)
  invisible(a)
}

How much crude oil does the world have available? Inquiring minds want to know but it is surprisingly difficult to get a single number. We hear more about relative levels “less than at any point since XXXX” and changes “has dropped X million barrels in the last month” than estimates of the absolute level. The International Energy Agency (IEA) publishes a lot of information on inventory levels and changes, but consolidated global stock estimates are, as far as I can see, not publicly available.

However, the Joint Organizations Data Initiative (JODI) does publish monthly country level estimates, reported to it voluntarily by countries. JODI partners are the APEC Energy Working Group, Eurostat, GECF, IEA, the Latin American and Caribbean Energy Organization (OLACDE), OPEC and the UN Statistical Division. It’s a very reputable set of partners! And JODI’s data are available to download.

Only 38 countries have data for every month JODI collects data, but this includes the USA, Japan, UK, the key European countries, and a few other big holders like Saudi Arabia, South Korea, Taiwan and Nigeria. The key omissions include, unfortunately but unsurprisingly, China, Russia, and Venezuela.. The data are also a bit old—latest data relating to July 2026 as at the time of writing, 6 October 2026, so a lag of 2-3 months. But it’s much better than nothing.

The point of today’s post is to get me the chart below, showing the best picture we can of the stocks of crude oil of the countries that do report those stocks, compared to a rough estimate of 85 million barrels per day of global crude refinery throughput:

That chart is now part of my fuel crisis monitoring page which gets updated at least weekly.

The overall story is that we indeed have less crude oil in storage than at any time since these records began in the early 2000s; but it’s still around 20 days of refinery cover. And at the rate that it’s going down, it will take nearly two years to get down to just 10 days of refinery cover. Actually, I’d expect market panic to set in at around 15 days of cover, but I’m being conservative here because I’m interested in the point where the world really has run out of any buffer against supply disruptions.

OK, the rest of the blog is mostly about just how I got that data, made some choices about excluding some countries, and drew the final chart.

First, downloading the data. There’s a separate CSV for each year, going back to 2002. Here’s how I download all those and import them into R

library(tidyverse)
library(glue)
library(janitor)
library(countrycode)

# Set to FALSE if you don't want to download the latest 2026 data. Will need to
# update all this workflow once we get to 2027 data.
update_2026 <- TRUE

#-----------------Downloads----------------

if (update_2026) {
  dir.create("fuel-crisis", showWarnings = FALSE)
  # Download current year. This will get more complete month by month.
  df <- here("fuel-crisis/jodi-2026.csv")
  url <- "https://www.jodidata.org/_resources/files/downloads/oil-data/annual-csv/primary/primaryyear2026.csv"
  download.file(url, destfile = df)
}


# Download historical data; only needed to be done once
if (!file.exists(here("fuel-crisis/jodi-2025.csv"))) {
  for (y in 2002:2025) {
    df <- here(glue("fuel-crisis/jodi-{y}.csv"))
    url <- glue(
      "https://www.jodidata.org/_resources/files/downloads/oil-data/annual-csv/primary/{y}.csv"
    )
    download.file(url, destfile = df)
  }
}

#----------------Import data------------------

# Import all the data to R:
jodi_hist_l <- list()
for (i in 1:25) {
  df <- here(glue("fuel-crisis/jodi-{i+2001}.csv"))
  jodi_hist_l[[i]] <- read.csv(df) |> as_tibble()
}

jodi_hist <- bind_rows(jodi_hist_l) |>
  clean_names()

# see https://www.jodidata.org/_resources/files/downloads/oil-data/jodi-oil-wdb-item-names-ver2017.pdf
# for guide on what everything is
# CLOSTLV is closing stocks
# CRUDEOIL is just crude
# TOTCRUDE also includes NGL and refinery feedstocks, additives and other hydrocarbons and has a lot more missing data.

Now that I’ve got the data, I want to understand the missingness. At one point I contemplated imputing the missing values of various countries, and hence I made a big wide version of the data with the idea that I might model crude stocks based on movement in some of the other numerous variables in the data. I eventually abandoned this when I realised how little meaningful information there was on eg China and Russia to do this meaningfully. But I can still use the jodi_wide data set created below for most of my plots:

#-----------explore which countries missing---------------------
jodi_wide <- jodi_hist |>
  #  filter(energy_product == "CRUDEOIL") |>
  mutate(
    dimensions = paste(energy_product, flow_breakdown, unit_measure, sep = "|")
  ) |>
  mutate(obs_value = as.numeric(obs_value)) |>
  select(ref_area, time_period, dimensions, obs_value) |>
  drop_na() |>
  spread(dimensions, obs_value, fill = NA) |>
  rename(crude = `CRUDEOIL|CLOSTLV|KBBL`) |>
  complete(ref_area, time_period, fill = list(crude = NA)) |>
  mutate(
    country = countrycode(
      ref_area,
      origin = "iso2c",
      destination = "country.name.en"
    )
  ) |>
  mutate(
    ref_area = fct_relevel(ref_area, "US"),
    date = ym(time_period),
    date_n = as.numeric(date)
  )

last_date <- jodi_wide |>
  filter(!is.na(crude)) |>
  summarise(ld = max(date)) |>
  mutate(ld = format(ld, "%B %Y")) |>
  pull(ld)

the_caption <- glue(
  "Source: JODI-OIL World Database. Most recent data is for {last_date}."
)


jodi_wide |>
  ggplot(aes(x = date, y = crude / 1000, colour = ref_area)) +
  geom_line() +
  theme(legend.position = "none") +
  scale_y_continuous(label = comma) +
  labs(
    x = "",
    y = "Millions of barrels",
    caption = the_caption,
    title = "Stocks of crude oil",
    subtitle = "Coloured by country; legend not shown."
  )

Here’s that first exploratory plot, a line chart of all the countries. The big one here is the USA; I’ve not shown the legend because it would take up too much space.

Next I want to identify which countries have missing data, and in particular those that are potentially big contributors to world crude oil stocks but are just missing a few observations that perhaps I could impute. I start by making a country_sum summary dataset:

possible_obs <- length(unique(jodi_wide$time_period))

country_sum <- jodi_wide |>
  group_by(ref_area, country) |>
  summarise(
    crude_total = mean(crude, na.rm = TRUE),
    missing_obs = sum(is.na(crude))
  ) |>
  ungroup() |>
  # particular problem if a country has any missing observations and averages
  # 10 million barrels or more of stock:
  mutate(
    all_there = as.logical(missing_obs == 0),
    some_missing = as.logical(
      missing_obs > 0 & missing_obs < possible_obs & !is.na(crude_total)
    ),
    problem = as.logical(crude_total > 10000 & missing_obs > 0)
  ) |>
  arrange(desc(crude_total))

good_countries <- country_sum |>
  filter(all_there) |>
  pull(country)

good_countries

That gives me this list of the countries that have a complete set of observations (largest crude oil stocks listed first):

> good_countries
 [1] "United States"  "Japan"          "Saudi Arabia"  
 [4] "Germany"        "South Korea"    "Canada"        
 [7] "France"         "Turkey"         "Italy"         
[10] "Poland"         "Spain"          "United Kingdom"
[13] "Netherlands"    "Taiwan"         "Thailand"      
[16] "Norway"         "Nigeria"        "Australia"     
[19] "Sweden"         "Finland"        "Czechia"       
[22] "Hungary"        "Portugal"       "Austria"       
[25] "Slovakia"       "Belgium"        "Denmark"       
[28] "Chile"          "New Zealand"    "Azerbaijan"    
[31] "Brunei"         "Ireland"        "Switzerland"   
[34] "Estonia"        "Iceland"        "Luxembourg"    
[37] "Latvia"         "Slovenia"   

It also lets me find some countries that have no observations at all:

# 20 countries that never have any observations for Crude:
never_crude <- filter(country_sum, is.na(problem))
print(never_crude$country)

> print(never_crude$country)
 [1] "Albania"             "Bangladesh"         
 [3] "Bermuda"             "Belarus"            
 [5] "China"               "Egypt"              
 [7] "Georgia"             "Hong Kong SAR China"
 [9] "Moldova"             "Malta"              
[11] "Malaysia"            "Niger"              
[13] "Nepal"               "Sudan"              
[15] "Singapore"           "Syria"              
[17] "Eswatini"            "Tajikistan"         
[19] "Vietnam"             "Yemen" 

Discovering that there were 20 of these, including China and Singapore, is what made me decide there was no point trying to impute the missing values. But I did a bit more exploration of which countries had some observations but not a complete set. Firstly, a scatter plot of all such countries:

…and then a faceted time series of the 9 biggest problematic countries.

I judged that maybe I could impute values for India, Mexico , Venezuela and Brazil based on their levels going up or down on average at about the world rate but decided this just wasn’t going to be worth it. And clearly Russia dwarfs them all and is just not worth trying, given the special circumstances there.

Here’s the code for those two plots:

country_sum |>
  ggplot(aes(
    x = crude_total,
    y = missing_obs,
    label = ref_area,
    colour = problem
  )) +
  geom_text() +
  scale_x_log10(label = comma) +
  labs(
    x = "Average crude stocks (thousands of barrels)",
    y = "Number of months missing an observation",
    caption = the_caption,
    title = "Size of crude stocks by number of missing observations"
  ) +
  theme(legend.position = "none")

# The nine biggest problem countries in terms of partly missing data
jodi_wide |>
  #  filter(ref_area %in% filter(country_sum, some_missing)$ref_area) |>
  filter(ref_area %in% filter(country_sum, problem)$ref_area) |>
  mutate(country = fct_reorder(country, crude)) |>
  ggplot(aes(x = date, y = crude / 1000)) +
  facet_wrap(~country, scales = "free_y") +
  geom_line(colour = "steelblue") +
  expand_limits(y = 0) +
  scale_y_continuous(label = comma) +
  labs(
    x = "",
    y = "Millions of barrels",
    title = "Crude oil stocks of countries missing at least one data point",
    caption = the_caption
  )

One last exploratory plot to try to see the impact of including or dropping those countries with partial data, here’s an area chart:

Thinking of the Russia situation in particular, I decided best just to exclude those partially missing countries altogether.

Here’s the code for that chart (and which also makes the crude_stocks data object I will shortly be using for my real plot)

#----------------Summarise and draw chart------------

# Crude oil stocks (excludes NGL etc because often missing data):
crude_stocks <- jodi_wide |>
  left_join(select(country_sum, ref_area, some_missing), by = "ref_area") |>
  group_by(date, some_missing) |>
  summarise(
    total_crude_mbbl = sum(crude, na.rm = TRUE) / 1000,
    reporting_countries = length(unique(ref_area))
  ) |>
  ungroup()

crude_stocks |>
  ggplot(aes(x = date, y = total_crude_mbbl, fill = some_missing)) +
  geom_area() +
  labs(
    caption = the_caption,
    y = "Total crude oil (millions of barrels)",
    x = "",
    title = "Total crude oil by missingness status of countries",
    fill = "Countries missing any data:"
  )

Finally, the code to draw the actual presentation plot. A key magic number here is the world refinery usage of crude oil at 85 million barrels per day; this is the approximate “crude runs” value from a recent IEA Oil Market Report. Note that this is less than the commonly cited values of around 105 million barrels per day, which includes NGLs, biofuels, and other liquid fuels. The 85m per day value is the appropriate one for me to use when looking at the CRUDEOIL series in my original JODI data; noting that they also have a TOTCRUDE series that would relate to the 105m per day figure, but which has a lot more missing data and is arguably not really what I’m most interested in anyway.

# note the commonly used figure of 105 is for TOTCRUDE, not just CRUDEOIl (which is more like 85)
world_use_per_day <- 85

# Summary for calculation of how long until only 10 days of cover:
war_summary <- crude_stocks |>
  # restrict it to countries that have observations for all months
  filter(!some_missing) |>
  summarise(
    prewar = total_crude_mbbl[date == "2026-03-01"],
    latest = total_crude_mbbl[date == max(date)],
    days = max(date) - as.Date("2026-03-01")
  ) |>
  mutate(
    rate = (prewar - latest) / as.numeric(days),
    ten_days_left = (latest - world_use_per_day * 10) / rate
  )

# draw chart:
crude_stocks |>
  filter(!some_missing) |>
  ggplot(aes(x = date, y = total_crude_mbbl)) +
  geom_line(colour = "blue") +
  geom_hline(yintercept = world_use_per_day * 10, colour = "darkred") + #expand_limits(y = world_use_per_day * 10) +
  expand_limits(y = 0) +
  scale_y_continuous(
    label = comma_format(suffix = "m"),
    sec.axis = sec_axis(
      ~ . / world_use_per_day,
      name = glue(
        "Days of crude inventory cover,\nat {world_use_per_day} million barrels per day"
      )
    )
  ) +
  labs(
    x = glue(
      "{length(good_countries)} countries in total have data for all months in this period."
    ),
    title = "Crude oil stocks of consistently reporting countries worldwide",
    y = "Millions of barrels",
    subtitle = glue(
      "Excluding countries with any missing data (eg China, Russia, India, Venezuela).
If the rate of decline since March 2026 continued, crude inventories would fall to 10 days of cover in {round(war_summary$ten_days_left / 7)} weeks."
    ),
    caption = the_caption
  )

That gets me this graphic, which is the one I will keep up to date from now on:

That’s all for today. Take care out there, and buy yourself an electric vehicle and rooftop solar if you can!

This post presents a function mlReserve() for machine-learning loss reserving on ‘real-world’ triangles, based on a fork of R’s ‘ChainLadder’ package (https://github.com/thierrymoudiki/ChainLadder/tree/master). It includes two experiments, each scored against a known true reserve: (A) Two triangles from an individual-claims simulator (Wang & Wüthrich, https://github.com/actuarial-data-science/PackageIndividualClaimsSimulator), embedded below: chain ladder vs mlReserve, in detail; (B) A mini-benchmark on freshly simulated triangles from a simpler aggregate generator, with tunable distortions.

# =============================================================================
# mlReserve: machine-learning loss reserving on "real-world" triangles
#
# Colab: Runtime > Change runtime type > R
#
# Two experiments, each scored against a known true reserve:
#   A. Two triangles from an individual-claims simulator (Wang & Wüthrich),
#      embedded below: chain ladder vs mlReserve, in detail.
#   B. A mini-benchmark on freshly simulated triangles from a simpler aggregate
#      generator (defined in step 4), with tunable distortions.
# =============================================================================

# ---- 0. Setup -------------------------------------------------------------------

options(Ncpus = 2, repos = "https://cloud.r-project.org",
        repr.plot.width = 12, repr.plot.height = 5.5)
for (pkg in c("remotes", "ggplot2", "ranger", "e1071"))
  if (!requireNamespace(pkg, quietly = TRUE)) install.packages(pkg)
if (!requireNamespace("ChainLadder", quietly = TRUE) ||
    !"mlReserve" %in% getNamespaceExports("ChainLadder"))
  remotes::install_github("thierrymoudiki/ChainLadder", upgrade = "never")

suppressPackageStartupMessages({
  library(ChainLadder)
  library(ggplot2)
})

n      <- 10                                   # triangle size (accident years)
models <- c("Chain ladder", "mlReserve: SVM", "mlReserve: random forest")
ink    <- "#0b0b0b"
colours <- c("Chain ladder"             = "#eb6834",
             "mlReserve: SVM"           = "#2a78d6",
             "mlReserve: random forest" = "#1baf7a",
             "Truth"                    = ink)
theme_set(theme_minimal(base_size = 13) +
            theme(legend.position = "top", legend.title = element_blank(),
                  panel.grid.minor = element_blank(), plot.title.position = "plot",
                  strip.text = element_text(face = "bold")))


# ---- 1. Helpers -----------------------------------------------------------------

# A full n x n square of incremental payments -> observed (cumulative) upper
# triangle + true outstanding reserve per accident year (sum of the lower triangle)
make_case <- function(full) {
  dimnames(full) <- list(origin = 1:n, dev = 1:n)
  future <- row(full) + col(full) > n + 1
  upper  <- full; upper[future] <- NA
  list(full      = full,
       upper     = incr2cum(as.triangle(upper)),
       true_ibnr = rowSums(full * future))
}

# Fit the three models; `uncertainty = FALSE` gives point estimates only
fit_models <- function(triangle, uncertainty = TRUE, nsim = 300) {
  mse <- if (uncertainty) "bootstrap" else "none"
  fits <- list(
    MackChainLadder(triangle, est.sigma = "Mack"),
    mlReserve(triangle, e1071::svm, features = "numeric", transform = "asinh",
              mse.method = mse, nsim = nsim, seed = 1),
    mlReserve(triangle, ranger::ranger, num.trees = 300, num.threads = 1,
              features = "both", transform = "asinh",
              mse.method = mse, nsim = nsim, seed = 1))
  setNames(fits, models)
}

# Reserve (IBNR) and standard error by accident year, whatever the model
reserve_by_origin <- function(fit) {
  if (inherits(fit, "MackChainLadder")) {
    s <- summary(fit)$ByOrigin
    data.frame(origin = 1:n, IBNR = s$IBNR, SE = s$Mack.S.E)
  } else {
    s <- fit$summary[rownames(fit$summary) != "total", ]
    data.frame(origin = as.integer(rownames(s)), IBNR = s$IBNR, SE = s$S.E)
  }
}
total_reserve <- function(fit) sum(reserve_by_origin(fit)$IBNR)


# ---- 2. Data for experiment A ---------------------------------------------------

# Incremental paid (thousands), aggregated from the Wang & Wuthrich individual
# claims simulator (github.com/actuarial-data-science/PackageIndividualClaimsSimulator):
#   "Speed-up"     : settlement accelerates for recent accident years (delays -45%)
#   "All combined" : speed-up + calendar inflation shock (2% -> 15% a year from
#                    calendar month 78) + a few large, slowly paid claims
squares <- list(
  "Speed-up" = c(4172.1, 16153.4, 14800.2, 7561.1, 3275, 1856.3, 726.7, 410.7, 35.8, 140, 4983.2, 18064.1, 13574.5, 8350.7, 3954.9, 1499, 510.3, 9.5, 0, 11.3, 4698.1, 19428.8, 14864, 7343, 3710.5, 1378.6, 384.7, 92.2, 6.8, 27.1, 5467.1, 18481, 14847.6, 8089.7, 2794.1, 1379.2, 549.5, 152.2, 6.8, 0, 5985.2, 21516.6, 15980.6, 8353, 2954.2, 1655.5, 636.4, 191.3, 211.3, 28.9, 8274.1, 22014.4, 15981.4, 7425.8, 2362.7, 731, 242.1, 41.6, 52.8, 0, 8148.4, 24712.9, 16329.1, 6587.5, 2804.9, 753.6, 49.5, 0, 0, 0, 9974.2, 27573.4, 14523, 5141.6, 1123.8, 350, 73.3, 0, 0, 0, 10698.8, 29877.8, 14303.3, 3200.3, 978.3, 239.1, 0, 0, 0, 0, 13669.1, 30541.7, 10984, 1837.5, 157.3, 0, 0, 0, 0, 0),
  "All combined" = c(4236.8, 16766.9, 15671.5, 8229, 3582.8, 2143.9, 838.8, 548.6, 57.6, 261.6, 5195.7, 19008.8, 15926.9, 9857.6, 4643.8, 1746.5, 712.4, 26.4, 0, 22.5, 4956.8, 20861.6, 16283.4, 8704.2, 4321.7, 1876.2, 666.5, 161.8, 14.7, 61.4, 5908.2, 20336.4, 16580.1, 9331.1, 3673.8, 2276.2, 1205.8, 332.2, 15.4, 0, 6586.1, 24045.6, 18499, 11218.3, 4576, 3174.8, 1307, 460.6, 665, 96.3, 9292, 26084.3, 21055.7, 11396.7, 4203.8, 1531.8, 585.4, 115, 163.7, 0, 9702.9, 33010, 25202.7, 12133, 6404.4, 1956.8, 385.9, 0, 0, 0, 13752.5, 42686.8, 25857.7, 10698.4, 2706.5, 3213, 4232.5, 0, 0, 0, 17094.3, 53962.4, 30132.7, 8527.6, 3087.9, 746.2, 0, 0, 0, 0, 25375.8, 63473.7, 26195.3, 5737.2, 507, 674.2, 0, 0, 0, 0))
cases <- lapply(squares, function(v) make_case(matrix(v, n, n, byrow = TRUE)))
combo <- cases[["All combined"]]

cat("Observed cumulative paid, 'All combined' (thousands):\n")
print(round(combo$upper))
fits <- suppressWarnings(fit_models(combo$upper))
cat("\nmlReserve with a random forest; true total reserve =",
    format(round(sum(combo$true_ibnr)), big.mark = ","), "\n")
print(fits[["mlReserve: random forest"]])


# ---- 3. Experiment A: one triangle in detail ------------------------------------

# 3a. Reserve by accident year vs the truth
by_origin <- do.call(rbind, lapply(models, function(m)
  data.frame(model = m, reserve_by_origin(fits[[m]]))))
by_origin <- subset(by_origin, origin >= 2)             # AY 1 is fully developed
by_origin$model <- factor(by_origin$model, levels = models)
truth_points <- data.frame(origin = 2:n, IBNR = combo$true_ibnr[2:n])

print(
  ggplot(by_origin, aes(origin, IBNR / 1000, colour = model)) +
    geom_linerange(aes(ymin = (IBNR - SE) / 1000, ymax = (IBNR + SE) / 1000),
                   position = position_dodge(0.6), linewidth = 0.9) +
    geom_point(position = position_dodge(0.6), size = 2.6) +
    geom_point(data = truth_points, aes(origin, IBNR / 1000), inherit.aes = FALSE,
               shape = 23, size = 3.4, fill = ink, colour = "white") +
    scale_colour_manual(values = colours) +
    scale_x_continuous(breaks = 2:n) +
    labs(title = "Reserve by accident year, +/- one standard error (black diamonds: truth)",
         x = "Accident year", y = "Reserve (millions)"))

# 3b. Predictive distribution of the total reserve. Each model uses its own
#     resampling scheme: ODP bootstrap (BootChainLadder) for chain ladder,
#     residual bootstrap (mlReserve) for the learners.
set.seed(1)
boot <- BootChainLadder(combo$upper, R = 2000, process.distr = "od.pois")
draws <- rbind(
  data.frame(model = models[1], total = boot$IBNR.Totals),
  data.frame(model = models[2], total = rowSums(fits[[2]]$sims.reserve.pred)),
  data.frame(model = models[3], total = rowSums(fits[[3]]$sims.reserve.pred)))
draws$model <- factor(draws$model, levels = models)

print(
  ggplot(draws, aes(total / 1000, colour = model)) +
    geom_density(linewidth = 1, adjust = 1.2) +
    geom_vline(xintercept = sum(combo$true_ibnr) / 1000, colour = ink, linewidth = 0.8) +
    annotate("text", x = sum(combo$true_ibnr) / 1000, y = Inf, label = " true reserve",
             hjust = 0, vjust = 1.5, colour = ink) +
    scale_colour_manual(values = colours) +
    labs(title = "Model-specific predictive distributions of the total reserve",
         x = "Total reserve (millions)", y = "Density"))

cat("\nPredictive quantiles of the total reserve (thousands):\n")
quants <- t(sapply(split(draws$total, draws$model), quantile, c(0.025, 0.5, 0.975)))
print(noquote(formatC(round(quants), format = "d", big.mark = ",")))

# 3c. Projected cumulative payments vs the truth, accident years 7-10
cumulative_paths <- function(case, label) {
  cl <- MackChainLadder(case$upper, est.sigma = "Mack")
  rf <- mlReserve(case$upper, ranger::ranger, num.trees = 300, num.threads = 1,
                  features = "both", transform = "asinh", mse.method = "none", seed = 1)
  projections <- list("Truth"                    = t(apply(case$full, 1, cumsum)),
                      "Chain ladder"             = unclass(cl$FullTriangle),
                      "mlReserve: random forest" = unclass(rf$FullTriangle))
  out <- expand.grid(dev = 1:n, origin = 7:n, model = names(projections),
                     stringsAsFactors = FALSE)
  out$value <- mapply(function(o, d, m) projections[[m]][o, d],
                      out$origin, out$dev, out$model)
  out$panel <- paste0(label, " · AY ", out$origin)
  subset(out, dev >= n + 1 - origin)      # from the last observed diagonal onwards
}
paths <- do.call(rbind, Map(cumulative_paths, cases, names(cases)))
paths$model <- factor(paths$model, levels = c("Truth", models[c(1, 3)]))
paths$panel <- factor(paths$panel, levels = unique(paths$panel))

options(repr.plot.height = 7)
print(
  ggplot(paths, aes(dev, value / 1000, colour = model, linetype = model)) +
    geom_line(linewidth = 1) +
    facet_wrap(~ panel, nrow = 2, scales = "free_y") +
    scale_colour_manual(values = colours) +
    scale_linetype_manual(values = c("solid", "22", "solid")) +
    scale_x_continuous(breaks = seq(2, n, 2)) +
    labs(title = "Projected cumulative payments from the last observed diagonal",
         x = "Development year", y = "Cumulative paid (millions)"))
options(repr.plot.height = 5.5)


# ---- 4. Experiment B: mini-benchmark on fresh triangles -------------------------

# Aggregate generator: gamma payment pattern with over-dispersed Poisson noise.
#   speedup   : mean payment delay shrinks by up to this share for recent AYs
#   inflation : calendar-year log-growth rate after calendar period 7 (2% before)
simulate_case <- function(seed, speedup = 0.3, inflation = 0.08, phi = 5) {
  set.seed(seed)
  ultimate <- 5e4 * (1 + 0.03 * (0:(n - 1))) * exp(rnorm(n, 0, 0.05))
  delay    <- 3 * (1 - speedup * pmax(0, (1:n - 4) / (n - 4)))
  pattern  <- t(sapply(delay, function(m) diff(pgamma(0:n, shape = 2, scale = m / 2))))
  calendar <- outer(1:n, 1:n, "+") - 1
  infl     <- exp(0.02 * pmin(calendar, 7) + inflation * pmax(calendar - 7, 0))
  mean_incr <- ultimate * pattern * infl
  make_case(phi * matrix(rpois(n * n, mean_incr / phi), n))
}

n_rep <- 20                                   # triangles per scenario
scenarios <- list(
  # "Baseline (no distortion)" = c(speedup = 0,   inflation = 0.02),  # control
  "Settlement speed-up"      = c(speedup = 0.3, inflation = 0.02),
  "Speed-up + inflation"     = c(speedup = 0.3, inflation = 0.08))

bench <- do.call(rbind, lapply(names(scenarios), function(sc) {
  do.call(rbind, lapply(seq_len(n_rep), function(seed) {
    case <- simulate_case(seed, speedup   = scenarios[[sc]][["speedup"]],
                                inflation = scenarios[[sc]][["inflation"]])
    point_fits <- suppressWarnings(fit_models(case$upper, uncertainty = FALSE))
    data.frame(scenario = sc, seed = seed, model = models,
               error = sapply(point_fits, total_reserve) - sum(case$true_ibnr))
  }))
}))

score <- aggregate(error ~ scenario + model, data = bench, FUN = function(e)
  c(bias = mean(e), RMSE = sqrt(mean(e^2)), MAE = mean(abs(e))))
score <- do.call(data.frame, score)
names(score) <- c("scenario", "model", "bias", "RMSE", "MAE")
score <- score[order(score$scenario, score$RMSE), ]

cat("\nError of the total reserve over", n_rep, "triangles per scenario (thousands):\n")
fmt <- function(x) formatC(round(x), format = "d", big.mark = ",")
print(transform(score, bias = fmt(bias), RMSE = fmt(RMSE), MAE = fmt(MAE)),
      row.names = FALSE)

long <- rbind(data.frame(score[c("scenario", "model")], metric = "RMSE", value = score$RMSE),
              data.frame(score[c("scenario", "model")], metric = "MAE",  value = score$MAE))
long$model <- factor(long$model, levels = rev(models))
print(
  ggplot(long, aes(value / 1000, model, fill = model)) +
    geom_col(width = 0.6) +
    facet_grid(metric ~ scenario, scales = "free_x") +
    scale_fill_manual(values = colours, guide = "none") +
    labs(title = "Out-of-sample error of the total reserve (lower is better)",
         x = "Millions", y = NULL))

Observed cumulative paid, 'All combined' (thousands):
      dev
origin     1     2     3     4     5     6     7     8     9    10
    1   4237 21004 36675 44904 48487 50631 51470 52018 52076 52338
    2   5196 24204 40131 49989 54633 56379 57092 57118 57118    NA
    3   4957 25818 42102 50806 55128 57004 57670 57832    NA    NA
    4   5908 26245 42825 52156 55830 58106 59312    NA    NA    NA
    5   6586 30632 49131 60349 64925 68100    NA    NA    NA    NA
    6   9292 35376 56432 67829 72032    NA    NA    NA    NA    NA
    7   9703 42713 67916 80049    NA    NA    NA    NA    NA    NA
    8  13752 56439 82297    NA    NA    NA    NA    NA    NA    NA
    9  17094 71057    NA    NA    NA    NA    NA    NA    NA    NA
    10 25376    NA    NA    NA    NA    NA    NA    NA    NA    NA

mlReserve with a random forest; true total reserve = 174,050 
mlReserve (formula interface)

        Latest Dev.To.Date Ultimate   IBNR         S.E        CV
2      57118.1   0.9994925  57147.1     29    39.76128 1.3710786
3      57832.2   0.9982222  57935.2    103    91.80065 0.8912685
4      59311.6   0.9935674  59695.6    384   244.97640 0.6379594
5      68099.8   0.9732881  69968.8   1869   888.80702 0.4755522
6      72032.5   0.9286787  77564.5   5532  1605.55443 0.2902304
7      80048.6   0.8558364  93532.6  13484  2443.15573 0.1811892
8      82297.0   0.7087359 116118.0  33821  7158.48746 0.2116581
9      71056.7   0.5017785 141609.7  70553 13547.27988 0.1920156
10     25375.8   0.1881838 134845.8 109470 23912.17707 0.2184359
total 573172.3   0.7090055 808417.3 235245 37338.60673 0.1587222

Predictive quantiles of the total reserve (thousands):
                         2.5%    50%     97.5%  
Chain ladder             316,100 356,087 397,519
mlReserve: SVM           178,587 252,035 371,639
mlReserve: random forest 162,579 217,369 301,946

Error of the total reserve over 20 triangles per scenario (thousands):
             scenario                    model    bias    RMSE     MAE
  Settlement speed-up mlReserve: random forest  67,436  67,519  67,436
  Settlement speed-up             Chain ladder  94,022  94,097  94,022
  Settlement speed-up           mlReserve: SVM  93,888  94,318  93,888
 Speed-up + inflation mlReserve: random forest  58,834  58,975  58,834
 Speed-up + inflation           mlReserve: SVM  79,344  79,800  79,344
 Speed-up + inflation             Chain ladder 100,217 100,311 100,217

The experiment is not a universal proof that random forests are better. It is a demonstration that mlReserve() can be used to explore the performance of machine-learning models on loss reserving problems, and that it can outperform the classical chain ladder method in some scenarios.

Pipeline setup

{pipeflow} builds a pipeline by adding one R function per step with pip_add (see also the Get started article). In this post, we’ll work with the following toy pipeline:

• step load produces a sequence
• step clean doubles it
• step fit sums it, and
• step report formats the result

pip <- pip_new("my-pip") |>
    pip_add(
        step = "load",
        fun = \(n = 5) seq_len(n),
        tags = c("io", "daily")
    ) |>
    pip_add(
        step = "clean",
        fun = \(x = ~load) x * 2,
        tags = c("io", "core")
    ) |>
    pip_add(
        step = "fit",
        fun = \(x = ~clean) sum(x),
        tags = c("model", "core")
    ) |>
    pip_add(
        step = "report",
        fun = \(x = ~fit) paste("result:", x),
        tags = "report"
    )

Each step must have a unique name, and function parameters can refer to other steps — for example, x = ~load takes the output of the load step. We also assign tags to each step, which we will use to filter steps later.

Before we run it, let’s have a look:

pip
<pipeflow> my-pip (4 steps)
---------------------------
     step params depends state       tags
1:   load      n           new   io,daily
2:  clean      x    load   new    io,core
3:    fit      x   clean   new model,core
4: report      x     fit   new     report
---------------------------
<ready> last run: never

We see one row per step, showing its parameters (params), the steps it depends on (depends), its current state (all new before the first run), and our tags. Once we run the pipeline, …

pip_run(pip)
info [2026-10-03 09:36:44.793 UTC]: Starting run of pipeflow 'my-pip'
info [2026-10-03 09:36:44.794 UTC]: Step 1/4 load
info [2026-10-03 09:36:44.795 UTC]: Step 2/4 clean
info [2026-10-03 09:36:44.797 UTC]: Step 3/4 fit
info [2026-10-03 09:36:44.799 UTC]: Step 4/4 report
info [2026-10-03 09:36:44.803 UTC]: Finished run of pipeflow 'my-pip'
pip
<pipeflow> my-pip (4 steps)
---------------------------
     step params depends state            out       tags
1:   load      n          done      1,2,3,4,5   io,daily
2:  clean      x    load  done  2, 4, 6, 8,10    io,core
3:    fit      x   clean  done             30 model,core
4: report      x     fit  done     result: 30     report
---------------------------
<ready> last run: 2026-10-03 11:36:44

… there is a new out column with the results, which can also be accessed directly:

pip[["report", "out"]]
[1] "result: 30"

Pipeline views in action

Real pipelines get long, and often you only care about a subset of steps: a topic, a stage, or the steps that produce the desired outputs. Views (new in v0.4.0) let you work on such a subset without copying anything. They reference the underlying pipeline, so every operation applied to it — running it, updating parameters, collecting output, … — writes through to the original pipeline, restricted to the steps the view covers.

Creating views

pip_view() returns a view containing only the steps that match the given filters:

pip_view(pip, tags = "core")
<pipeflow_view> my-pip view (2 of 4 steps)
------------------------------------------
    step params depends state            out       tags
1: clean      x    load  done  2, 4, 6, 8,10    io,core
2:   fit      x   clean  done             30 model,core
------------------------------------------
<ready> last run: 2026-10-03 11:36:44
pip_view(pip, tags = "core", state = "done")
<pipeflow_view> my-pip view (2 of 4 steps)
------------------------------------------
    step params depends state            out       tags
1: clean      x    load  done  2, 4, 6, 8,10    io,core
2:   fit      x   clean  done             30 model,core
------------------------------------------
<ready> last run: 2026-10-03 11:36:44
pip_view(pip, step = c("clean", "fit"))
<pipeflow_view> my-pip view (2 of 4 steps)
------------------------------------------
    step params depends state            out       tags
1: clean      x    load  done  2, 4, 6, 8,10    io,core
2:   fit      x   clean  done             30 model,core
------------------------------------------
<ready> last run: 2026-10-03 11:36:44

By default steps must match all filters (logical AND), while the values within one filter are alternatives (OR). join = "union" keeps steps that match any filter, and fixed = FALSE treats filter values as regular expressions:

pip_view(pip, tags = "report", step = "clean", join = "union")
<pipeflow_view> my-pip view (2 of 4 steps)
------------------------------------------
     step params depends state            out    tags
1:  clean      x    load  done  2, 4, 6, 8,10 io,core
2: report      x     fit  done     result: 30  report
------------------------------------------
<ready> last run: 2026-10-03 11:36:44
pip_view(pip, step = "^f", fixed = FALSE)
<pipeflow_view> my-pip view (1 of 4 steps)
------------------------------------------
   step params depends state out       tags
1:  fit      x   clean  done  30 model,core
------------------------------------------
<ready> last run: 2026-10-03 11:36:44

Selecting steps with [

The extract operator (also new in v0.4.0) provides a data.table-like way of selecting steps and returns a view by default:

pip[c("load", "fit")]
<pipeflow_view> my-pip view (2 of 4 steps)
------------------------------------------
   step params depends state       out       tags
1: load      n          done 1,2,3,4,5   io,daily
2:  fit      x   clean  done        30 model,core
------------------------------------------
<ready> last run: 2026-10-03 11:36:44

Boolean filters are evaluated in the context of the step table.

pip[tags %like% "core"]
<pipeflow_view> my-pip view (2 of 4 steps)
------------------------------------------
    step params depends state            out       tags
1: clean      x    load  done  2, 4, 6, 8,10    io,core
2:   fit      x   clean  done             30 model,core
------------------------------------------
<ready> last run: 2026-10-03 11:36:44
pip[step %in% c("clean", "fit") & state == "done"]
<pipeflow_view> my-pip view (2 of 4 steps)
------------------------------------------
    step params depends state            out       tags
1: clean      x    load  done  2, 4, 6, 8,10    io,core
2:   fit      x   clean  done             30 model,core
------------------------------------------
<ready> last run: 2026-10-03 11:36:44

Running views

Running a view executes the steps it covers together with any upstream dependencies that are not up to date. The run log marks the view’s own steps as [view] and the dependencies pulled in as [upstream]:

pip_reset(pip)  # reset the pipeline to "new" state for demonstration
pip_run(pip_view(pip, step = "report"))
info [2026-10-03 09:36:44.883 UTC]: Starting run of pipeflow 'my-pip view'
info [2026-10-03 09:36:44.883 UTC]: Step 1/4 [upstream] load
info [2026-10-03 09:36:44.883 UTC]: Step 2/4 [upstream] clean
info [2026-10-03 09:36:44.884 UTC]: Step 3/4 [upstream] fit
info [2026-10-03 09:36:44.885 UTC]: Step 4/4 [view] report
info [2026-10-03 09:36:44.886 UTC]: Finished run of pipeflow 'my-pip view'

Afterwards, the original pipeline is up to date for the covered steps. Let’s collect all results related to the “core” steps:

pip[tags %like% "core"] |> pip_collect()
$clean
[1]  2  4  6  8 10

$fit
[1] 30

For a complete guide on composing views, see the Pipeline views vignette.

More new and updated features

Views are only one of the additions. Each of the following has its own article on the documentation site:

• Functional and method-style API — the pipeable pip_*() functions, with every function also available as a method (p$run(), p$view(), …): Get started.
• Extract and replace operators — data.table-style [ and [[ for reading and editing steps: reference.
• Collect and group output — pip_collect() with grouping by tags or dependencies: Collect and group output.
• Combining and modifying pipelines — combine several pipelines with rbind(), or replace, rename, remove, and reset steps: Combine pipelines · Modify existing pipelines.
• Split, map, and reduce — native exec = "split", "auto", and "reduce" modes: Split, map, and reduce.
• Nested pipelines — use a pipeline as a reusable building block inside a step: Nested pipelines.
• Self-modifying pipelines — .self, run-time restructuring, and restart(): Self-modifying pipelines.

Under the hood, significant performance gains have been made by implementing the dependency graph in C++ (since v0.3.0) and by improving the step and parameter bookkeeping in the underlying {data.table}. Finally, by dropping lgr and jsonlite, external package dependencies have been reduced to {data.table} and Rcpp.

{pipeflow} vs {targets}

{targets} remains R’s de-facto standard for heavy-duty, reproducible pipelines: it persists results to disk, records provenance, and scales to distributed infrastructure via crew. {pipeflow} does not try to compete on that home turf but targets a different niche — the interactive session, often as the backend of a Shiny app — where a pipeline is built and modified on the fly and has to respond while a user waits. Since responsiveness is key in that use case, {pipeflow} has been optimized for low latency.

The pipeflow vs targets article compares the in-session latency of the two packages head to head. The figure below gives a flavour of the results:

It shows the scenario where one parameter has changed¹ and the pipeline is rerun, with each step doing 5 ms of work. The x-axis represents the number of steps in the pipeline (its size), and the y-axis shows the time required to re-run it.²

Across the tested sizes, {pipeflow} is consistently the faster of the two — see the article for the full benchmark and methodology.

Overall, the two are best seen as complementary: reach for {targets} for large batch jobs that must be reproducible end to end, and for {pipeflow} when the pipeline is part of an interactive application.

Wrapping up

{pipeflow} 0.4.0 is on CRAN:

install.packages("pipeflow")

The full changelog lists every change, and the documentation includes the articles linked above. If you build Shiny apps or spend a lot of time exploring parameters interactively, give it a try — feedback and issues are welcome on GitHub.

I was doing some analysis work with R today (don’t get to do that very often anymore) and R wasn’t playing with VS Code nicely.

Turns out there has recently been a major update to the VS Code R editor extension.

So if you are having issues, check out the Extensions page for setup instructions.

The North American Bioconductor Conference 2026 (([BioC2026])(https://bioc2026.bioconductor.org/)) took place from August 10-12, 2026, at the Fred Hutchinson Cancer Center in Seattle, Washington. The conference brought together researchers, developers, educators, and community members to explore current developments within and beyond the Bioconductor project. BioC2026 welcomed 108 attendees, with 70 participating in person and 38 joining virtually. Across the three days, the programme featured three keynote speakers, 20 contributed speakers, and 17 sessions. The programme highlighted the breadth of biological research supported by Bioconductor, as well as the tools, infrastructure, and communities that enable open and reproducible computational biology. The conference was then followed by two community events on August 13-14: a Carpentries workshop on orchestrating large-scale single-cell analysis with Bioconductor and the Bioconductor North America 2026 Hackathon.

Programme overview

Keynotes

The conference featured three keynote speakers: – Jeff Leek speaking on “Building Open Infrastructure for Translational AI in a Cancer Center” – Ting Ye speaking on “From Association to Causation: Genetic-Anchored Causal Inference in Human Biology” – Michael Lawrence speaking on “S7 for Bioconductor”

(Left to Right: Jeff Leek, Michael Lawrence and Ting Ye giving their keynotes at BioC2026)

Short talks

Twenty contributed speakers shared their work across the conference programme. The short talks at BioC2026 reflected the breadth of research and development across the Bioconductor community. Topics ranged from community training and developer engagement to multi-omics exploration, reproducible workflows, and AI infrastructure. The programme also covered emerging areas of biological data analysis, including spatial transcriptomics, antibody repertoire analysis, rare variant testing, DNA methylation, structural variation, and Perturb-seq.

Poster sessions and package demo sessions

The poster session showcased a diverse range of research applications and software development projects from across the Bioconductor community. Topics included single-cell and spatial transcriptomics, multi-omics integration, microbiome analysis, proteomics, antibody repertoire analysis, functional enrichment, ancestry inference, genomic data processing, and AI-guided drug discovery. The posters provided an opportunity for attendees to discuss new methods and tools in a more informal setting, while connecting researchers, package developers, and community members working across different areas of computational biology. The sessions also complemented the package demonstrations, creating space for participants to explore new software and engage directly with the people developing and applying Bioconductor tools.

Feedback from the conference showed that these interactions were particularly valuable, with attendees highlighting conversations around posters and package demonstrations as some of the most useful networking opportunities.

Participants during the poster session.

Celebrating 25 years of Bioconductor

The celebrations marking Bioconductor’s 25th anniversary continued at BioC2026 following the celebrations that began at EuroBioC2026 earlier in the year. The milestone was marked with a specially branded anniversary cake, thanking participants and community members for being part of the Bioconductor community and contributing to its growth over the past 25 years.

Bioconductor at 25 celebrations.

Post-conference events

Orchestrating Large-Scale Single-Cell Analysis with Bioconductor

On August 13–14, participants took part in a Carpentries workshop focused on orchestrating large-scale single-cell analysis with Bioconductor led by Ted Laderas and Wes Wilson. The workshop was attended by 11 participants and provided a hands-on environment for participants to work through approaches for analysing large-scale single-cell datasets using Bioconductor tools and workflows. Feedback highlighted the value of the live coding, small-group format and complementary teaching approaches.

Participants during the scRNA-seq workshop.

Bioconductor North America 2026 Hackathon

The Bioconductor North America 2026 Hackathon took place on August 13–14, bringing community members together to collaborate on software and infrastructure projects that strengthen interoperability, software development, and the wider Bioconductor ecosystem. These projects included: AI harness tools for standardizing how users within organizations would interact with agents on shared projects and on shared institutional resources. Link: https://github.com/eisenefaust/harness_agent

Bioconductor tooling for including CLI executable components in Bioconductor packages. Link: https://blog.bioconductor.org/posts/2026-08-21-biocjobs/

An exploration on benchmarking LLM generated workflows, building out infrastructure and paradigms for how we can evaluate what LLMs return as work product to scientific questions. Link: https://github.com/Amisor/LLM_BioBenchmarking/tree/main

Participants during the Hackathon.

© 2026 Bioconductor. Content is published under Creative Commons CC-BY-4.0 License for the text and BSD 3-Clause License for any code. | R-Bloggers

You can read the original post in its original format on Rtask website by ThinkR here: My first attempt at porting a Shiny app failed, and not for a technical reason

Diary of a rewrite, episode 2 of 4.

Translating your code line by line is the surest way to get it all wrong…

In the previous episode we told the story of why we cancelled our own stack decision halfway through, and why we went back to R. If you have not read it, no worries: this one stands on its own.

Today we are talking about the port itself. And above all about the way I got it wrong the first time.

A quick bit of context. StaffPuzzle is our internal staffing tool: it says who works on what, half-day by half-day. Its Shiny version is 1,200 lines, reads the team’s Google calendars, and has been running for years. The goal was to rewrite it in htmxr, meaning HTML rendered server-side, hydrated by htmx, with plumber2 behind it.

Here is how I went about it.

I opened app_server.R, I listed every observeEvent, and I started turning them into routes. One route per observer. It is methodical, it is exhaustive… and it produces an application full of endpoints that serve no purpose.

A quick reminder: what is a route? In a classic web application, the browser asks for an address and the server sends back a response. A route is that pair: an address, and the function that builds what we send back. Where Shiny holds a continuous conversation with the browser for the whole visit, a route-based application answers requests that are independent of one another, with no memory of the previous one.

Concretely, here is the kind of translation I was doing. On the Shiny side, the observer that rebuilds the main plot whenever a filter moves (I have shortened and renamed it for readability, the original is bushier):

observeEvent(c(input$refresh, rv$slots, rv$data_ready), {
  req(rv$slots, rv$consultants, rv$start, rv$end, rv$sort_by)
  rv$filtered_slots <- rv$slots |>
    filter(member %in% rv$consultants) |>
    mark_focus(types = rv$slot_types) |>
    add_weekdays()
  rv$plot <- rv$filtered_slots |>
    team_availability_plot(
      start   = rv$start,
      end     = rv$end,
      sort_by = rv$sort_by,
      order   = rv$consultants
    )
})

And, a hundred lines further down in the same file, the output that displays it:

output$plot <- renderPlotly({
  rv$plot |> ggplotly(tooltip = "text")
})

On the routes side, it becomes an address and a function:

#* @get /plot
function(consultants, start, end, sort_by, slot_types) {
  slots |>
    filter(member %in% consultants) |>
    mark_focus(types = slot_types) |>
    add_weekdays() |>
    team_availability_plot(
      start   = start,
      end     = end,
      sort_by = sort_by,
      order   = consultants
    )
}

The browser will ask for /plot?consultants=arthur,vincent&start=2026-10-05&end=2026-11-30&sort_by=alphabetical, and the server will answer with the corresponding piece of page.

Notice in passing what has melted away on the Shiny side: I need an observer, a slot in reactiveValues to store the result, and a render* elsewhere in the file to display it. On the route side, all that is left is a function that takes arguments and returns something. That one translates well: it produces something visible, and that something can state its address.

The problem is all the others.

My mistake did not come from not knowing the destination tool. I was translating code instead of changing paradigm. And it took me several days to notice!

It is a bit like moving house. You do not pack by photographing each room so you can rebuild it identically in the new flat. You look at what you actually use, and you leave the rest at the car boot sale.

So here is what I wish someone had told me before I started.

Why inventorying the observers is the wrong entry point

The reason is structural. Shiny is driven by a reactive graph, and a good share of the observers have no visible output at all. They copy a value from one place to another, they disable a button, they keep two inputs in sync. In short, they maintain the consistency of an in-memory state.

But that state does not exist in htmxr. Neither does the work of those observers. They do not translate: they disappear.

I remember the sentence I blurted out in a meeting, on 30 June, right in the middle of the port:

I feel like saying “I’ve got this or that observeEvent, so I’ve got this or that endpoint” isn’t ok. Whereas I thought starting from the data to serve is a good entry point.

It has been my working rule ever since.

Start from the data being served

An htmxr application is a set of addressable blocks of data. Every block on screen has a URL that renders it. The controls that drive it (filters, dates, sorting) are the parameters of that URL. We are in a hypermedia frame, no longer in a dependency graph.

The method comes down to 5 moves, and the order matters:

1. Inventory the screens, screenshot by screenshot. Not what the code does: what the user sees.
2. Cut each screen into blocks of data: header, filters, grid, detail panel.
3. Give every block a URL. A block that cannot state its address is a badly cut block.
4. Turn the controls into parameters. ?start=2026-10-01&n_weeks=8, not observers.
5. For each interaction, ask yourself a single question: is this a mutation of persistent business state?

The fifth one is the one I kept missing. Clicking a filter, changing a sort, picking a period: those are GETs, even if “something changes” on screen. Nothing is modified server-side, we are simply asking for another view of the same data.

Tip: if your inventory produces a lot of POSTs, chances are some GETs have slipped in. Go back over them.

60 lines with nothing left to do

Here is the example that tipped my understanding.

My Shiny version devoted 60 lines to keeping the URL and the inputs in sync, in both directions. Good code, by the way: that application had shareable URLs, which is not that common in Shiny.

On the Shiny side, 60 lines. Here are the two ends of it, abridged and renamed:

## Direction 1: the URL arrives, we put the inputs in the right state
observeEvent(session$clientData$url_search, {
  params <- parseQueryString(session$clientData$url_search)
  url_params$consultants <- params$consultants
  url_params$slot_types  <- params$types
  url_params$services    <- params$services
  url_params$period      <- params$period
  if (!is.null(url_params$consultants)) {
    updateSelectInput(
      session  = session,
      inputId  = "consultants",
      selected = str_to_title(
        str_split_1(str_remove_all(url_params$consultants, " "), pattern = ",")
      )
    )
  }
  ## … and the same treatment for period, slot_types and services
})
## Direction 2: an input moves, we rewrite the URL
observeEvent(
  c(input$consultants, input$slot_types, input$services,
    input$period, input$sort_by),
  ignoreInit = TRUE, {
    consultants <- str_to_lower(glue_collapse(input$consultants, sep = ","))
    slot_types  <- str_to_lower(glue_collapse(input$slot_types, sep = ","))
    period      <- glue_collapse(input$period, sep = ",")
    url <- glue("?consultants={consultants}&types={slot_types}&period={period}")
    updateQueryString(session = session, url, mode = "replace")
  })

On the htmxr side, one argument:

hx_set(
  form,
  get      = ".",
  target   = "#trainings-view",
  trigger  = "change",
  push_url = "true"
)

Careful, this is not a matter of concision. Those 60 lines have not been shortened: they have nothing left to do. They were maintaining a mirror between two representations of the same state, the URL on one side, the inputs on the other. When the URL is the state, there is simply no mirror left to hold.

Back to our house move: these are not furniture you carry up to the new flat. There is no room for them.

Want a simple test to see where you stand? If your port produces code that synchronises two things, you have kept the original paradigm. A good port deletes whole categories of code, it does not rewrite them shorter.

What happens to each building block

Once the frame is set, the port becomes mechanical. Let us look at what becomes of the blocks you have in front of you in your server.R.

A module becomes a route and a rendering function. My biggest module was 267 lines: 4 selectInput, one dateRangeInput, 6 observers, and the URL synchronisation. All that is left is an HTML form and some parameters. The ns() disappears entirely, since it only existed to avoid identifier collisions between instances. With no input binding, there is nothing left to name.

input$x becomes an argument of your function. There is no input object any more. The value arrives in the request: as a URL parameter for a GET, in the body for a POST. It is the most disorienting change in the first few days, and the most liberating one afterwards. A view function receives its inputs as arguments, so it can be tested!

output$x <- render*() becomes a function that returns HTML. The output/render pair disappears. All that is left is a function and its return value, called by a route. A renderPlot becomes a URL that serves the SVG, a renderUI a URL that serves a fragment.

updateSelectInput() becomes: nothing. We no longer update a control. The server renders the <select> in the state it should be in, and htmx swaps the fragment that contains it. The whole update*Input() family evaporates with it.

invalidateLater() becomes an attribute. A periodic refresh is written trigger = "every 15s". The browser asks for the fragment again, the server recomputes it. Nothing to invalidate, since there is no state.

showNotification() and withProgress() become a fragment, or CSS. A message is a piece of HTML returned with the response. A waiting indicator is an element that htmx shows during the request and hides afterwards. A CSS class, not R code.

reactiveValues() and session$userData become the URL, or the database. In the URL if it describes what we are looking at, in the database if it has to survive closing the tab. And if nothing is left after that sorting? Then it was plumbing state, and it goes with the rest

Seeing it run

All of this stays abstract until you have it in front of you. So let us take Shiny’s canonical example, the one you get with shiny::runExample("01_hello"): the Old Faithful geyser histogram, with a slider for the number of bins.

• The application, written in htmxr. One slider, one plot, and that is all.
• The URL that serves the plot. Open it directly in your browser: you land on the SVG itself, the drawing spelled out. That is exactly what “give every block a URL” means.

Have fun changing bins=30 to bins=5, then to bins=50. The server recomputes and sends you back the drawing. No session, no state, no open conversation.

And the code, on both sides

You already know the Shiny version: it is the one from the tutorial. Here it is, stripped to the essentials.

ui <- page_sidebar(
  title = "Hello Shiny!",
  sidebar = sidebar(
    sliderInput("bins", "Number of bins:", min = 1, max = 50, value = 30)
  ),
  plotOutput("distPlot")
)
server <- function(input, output) {
  output$distPlot <- renderPlot({
    x    <- faithful$waiting
    bins <- seq(min(x), max(x), length.out = input$bins + 1)
    hist(x, breaks = bins, col = "darkgray", border = "white")
  })
}

And the htmxr version, the one running behind the link above. I have removed the Bootstrap classes, they teach you nothing.

#* @get /
#* @serializer html
function() {
  hx_page(
    hx_head(title = "Old Faithful Geyser Data"),
    hx_slider_input(
      id      = "bins",
      label   = "Number of bins:",
      value   = 30, min = 1, max = 50,
      get     = "plot",
      trigger = "input changed delay:300ms",
      target  = "#plot"
    ),
    tags$div(id = "plot") |>
      hx_set(get = "plot", trigger = "load", target = "#plot")
  )
}
#* @get /plot
#* @serializer none
function(query) {
  xmlSVG({
    x    <- faithful$waiting
    bins <- seq(min(x), max(x), length.out = as.numeric(query$bins %||% 30) + 1)
    hist(x, breaks = bins, col = "darkgray", border = "white")
  }, width = 7, height = 5)
}

Three things to look at, and they sum up the whole article.

The slider is the same object. sliderInput(id, label, min, max, value) becomes hx_slider_input(id, label, min, max, value). There is nothing to relearn.

The plot computation has not moved by a single line. The four lines of hist() are identical on both sides. That is what I had not understood when I started: it is not your business code you are porting, it is the plumbing around it.

What changes is the wiring, and it becomes visible. In Shiny, nothing says that moving the slider redraws the plot: the reactive graph works it out by itself from the presence of input$bins inside renderPlot. In htmxr, it is written on the slider itself: get = "plot" (go fetch this address), target = "#plot" (put the result there). You read the wiring instead of guessing it.

And server <- function(input, output) is gone. All that is left is two functions that take arguments and return something.

Try it yourself

No need to copy anything: the example ships with the package.

# install.packages("pak")
pak::pak("hyperverse-r/htmxr")
htmxr::hx_run_example("hello")

Your browser is waiting at http://127.0.0.1:8080. Move the slider: every move fires a request to /plot, and the server returns a new SVG. Open your browser’s network tab to watch them go by, it beats any explanation.

And if you call hx_run_example() with no argument, it lists the others: select-input, multi-select, delete-row, toast-notification, infinity-scroll, page-or-fragment, json-endpoint and reactive-values. That last one is worth a look if you are wondering what becomes of a reactiveValues() once the reactive graph is gone.

And since the plot is a URL, anyone can consume it. Here is the same application written in React, calling exactly the same R endpoint. The backend does not know who is asking, and it does not care.

This is the most underrated property of cutting things into addressable blocks: it does not lock you in. The day you want another frontend, or a mobile app, or simply for a colleague to pull your plot into their own tool, your R logic is already an API. It became one without you doing anything extra.

What costs, and what does not

Everything above is mechanical. On my port, only 3 patterns called for a real decision, and they are the ones that make an estimate vary.

Shiny pattern	Cost	Equivalent
selectInput, actionButton, sliderInput	low	hx_select_input(), hx_button(), hx_slider_input()
dateRangeInput	low	no dedicated constructor: two <input type="date">
URL ︎ inputs sync	low	the push_url argument
static ggplot	low	svglite + htmx swap
DT (display)	low	hx_table(), hx_table_rows()
long task + notifications	medium	external job + periodic polling
interactive plotly	high	plotly.js, or drop the interactivity
rhandsontable	high	JS library, or one form per cell

The editable table has no native equivalent. Either you embed a JavaScript library, or you rebuild it with HTML fields and one save per cell. The criterion is dimensional: under twenty rows or so HTML wins, at a thousand rows with multi-cell copy-paste the library becomes the right choice again.

The interactive plot deserves a prior question, and it is precisely the one I forgot to ask myself: is it really a plot? When it is in fact a data grid, a schedule, a calendar or a heatmap, porting it as a plot is a mistake. Mine was one because Shiny made tables hard, not because the table was the right object. A real plot, on the other hand, stays an SVG served by a URL, and that is an easy case.

An opaque internal dependency is not worked around, it is isolated. You definitely have one: the homegrown package that fetches the data from somewhere, that nobody has really read in years. The temptation is to start with it, to “understand before porting”. Do not. You will spend three days on it and you will have ported nothing.

Develop everything else against a demo dataset, hard-coded, and only plug the real package in at the very end, in one go. By then you will know exactly what you want from it: a function, some arguments, a table back. And if it is still incomprehensible, call it as it is without trying to rewrite it. That is not cowardice, it is what stops it from blocking the whole port.

The most useful side effect

Shiny encourages mixing business logic and reactive orchestration: they live in the same file, often in the same function. In my v1, the rule that decides whether a slot is highlighted was computed twice, with different rules, 15 lines apart. That was not carelessness on my part, it is what a file where business and plumbing cohabit produces.

The port forces the separation. A function called by a route can no longer lean on an ambient reactive context. What is left is ordinary R, which takes data and returns data or HTML. So it can be tested.

The tests-to-code ratio went from 0.26 to 0.77: 544 tests, no browser, a few milliseconds. With the interface having become a pure function from data to HTML, it falls into the same test harness as the rest. The Shiny equivalent needs shinytest2, a headless Chrome and snapshots. It is slow, it is fragile, and it breaks when a margin moves.

Two pieces of honesty are in order here. My new version does more than the old one (writing to calendars, draft mode, historised indicators), and part of the difference in volume comes from that, not from the paradigm. And nothing stops you from testing the business logic of a Shiny application: v1 did it over 300 lines. What changes is that the interface stops being the area you never test.

What you do not get back

Let us be clear: the reactive graph does render a real service when everything depends on everything. Think of a dashboard where 6 linked plots recompute together, with cross dependencies and expensive derived values. Rebuilding that with fragments and URL parameters is possible, but you are hand-writing what Shiny gives you for free.

Hence the only decision criterion that is worth anything, to my mind: does my application’s state fit in a URL?

For mine, yes: a period, some filters, a sort. The port was a good bet. For an exploration tool where the user builds a complex working session, the answer is no, and Shiny remains the right tool. That is not a polite concession, it is the same reasoning applied to a different problem.

Conclusion

I still write Shiny, for our clients and for ourselves. This port did not prove that one technology beats another.

What it did show me is that for years I had been solving a page-rendering problem with a tool designed to hold state in memory. And that the cost of that gap stayed invisible until I had paid it once.

So before you set off, ask yourself the URL question. It will save you the few days I lost!

Do you have a Shiny application in this situation? We read it and give you back a note that says whether the port is worth it, including when the answer is no. Write to us in two lines: the size of your application, and how long it has been running.

Next week, episode 3: what leaving the managed platform cost us, bill in hand.

Diary of a rewrite, four episodes:

1. We chose React to learn it. The AI wrote it for us.
2. My first attempt at porting a Shiny app failed, and not for a technical reason (you are here)
3. We left the managed platform. Here is the real bill. (6 October)
4. Five signs your Shiny app has outgrown Shiny (13 October)

The application cited is StaffPuzzle, ThinkR’s internal staffing tool, written with htmxr, alpiner and plumber2. The figures are measured on the repository’s two branches, not estimated. Disclosure of interest: I am the author of htmxr and alpiner.

This post is better presented on its original ThinkR website here: My first attempt at porting a Shiny app failed, and not for a technical reason

Three hundred twenty-nine new packages were submitted to CRAN in July. Here are my Top 40 picks in nineteen categories: Agriculture, Artificial Intelligence, Climate Studies, Computational Methods, Ecology, Epidemiology, Genomics, Machine Learning, Medical Statistics, Meta-Analysis, Networks, Process control, Programming, Public Transit, Statistics, Surveys, Time Series, Utilities, and Visualization.

Read it in: Español.

rOpenSci turns 15 and I’ve been a part of this community for a good chunk of that time. As I reflected on those years, a lot of great memories came to mind, and I realized some of them taught me some lessons that I would love to share.

Stretch your comfort zone

Early in 2018 I came across this post:

Apply to attend rOpenSci unconf 2018!
We’re organizing unconf18 to bring together scientists, developers, and open data enthusiasts from academia, industry, government, and non-profits to get together for a couple of days to hack on various projects and generally enrich our community.

In hindsight I was a good fit. I had recently transitioned from academia to my first job as a research software engineer. And yet, did I think I was a solid candidate? Oh no, I had lots of doubts.

But I know the greatest things often happen a little stretch outside your comfort zone. So I applied anyway, got invited, and spent two days collaborating with an incredible bunch of people. Some of them were already rockstars, many others eventually became influential in their own niche, and all of them had this great attitude that at the time I lacked words to describe. I now do: they were respectful and kind.

Be respectful and kind

In 2022 rOpenSci launched the first cohort of the rOpenSci Champions Program. By then I already had a few years of experience as an associate editor of rOpenSci software peer-review, and Yani invited me to talk to our champions about the process.

As I prepared that talk and discussed some ideas with Yani, she made me realize how much the rOpenSci community cares about the way we communicate with one another¹. So much so that our guides for reviewers and editors open with this message:

rOpenSci’s community is our best asset. We aim for reviews to be open, non-adversarial, and focused on improving software quality. Be respectful and kind! See our reviewers’ guide and code of conduct for more.

Coming from such a deeply technical guide, that opening might surprise you. But it makes sense to me now; broken code is way easier to fix than broken human relationships.

Contribute in your own way

At rOpenSci you can contribute in so many ways. We have a community contributing guide but here’s a list of some of my own contributions:

• Attend events.
• Help to welcome and onboard new members.
• Ask or answer questions on Slack, or share or discuss ideas or jobs.
• Improve our documentation, e.g. fix a typo.
• Write a blog post.
• Translate or review some work in your native language.
• Pilot a new process.
• Author, review, or edit an R package.
• Lead a workshop.
• Build a tool to enhance some process.
• Give and get support, e.g. mentor, review grant applications or talks, advice or recommendations for jobs, help unavailable or overwhelmed people.

As you can see, most of my contributions did not involve any code, and only a few required an invitation. Some of them are the kind of thing you could add to your CV, and the most important ones only belong in your heart.

Attract great people

Great communities are made of great people, and you can play an active role in growing it in the direction you want. If you know someone that would fit in and enjoy the rOpenSci community, you can try to find a way to bring them in. For example, I’ve encouraged and helped people to submit and review packages. Also as an editor I have the occasional privilege to nominate other editors. Using this super-power I’ve attracted two editors and I now get to enjoy interacting with them quite regularly.

Keep the human in your human connections

AI is changing the way we contribute to open source software. But we don’t know much about its effect on the open communities around that software. Earlier this year rOpenSci published a preliminary set of policies, and you may want to watch rOpenSci’s newsletter and this Openscapes community call:

Open Communities in the Age of AI
{Open communities} help us connect on a human level around science and data and the things that we are passionate about. They change careers, they change lives. (…) But now, people are not necessarily finding their communities and getting the benefit of those communities in the way that they used to, because they’re more easily able to get help and answers by using AI. (…) These connections and communities are more important than ever, as we need to find ways to support each other and ourselves as AI rapidly changes the landscape right under our feet.

Personally, over the past year I’ve used AI heavily and learned several lessons. The most relevant here is that I don’t want AI to impersonate me or the other person in a human-to-human connection. The humans of rOpenSci are real and wonderful people. Here are some that I’ve recently met in person:

Thanks for joining me in celebrating rOpenSci. And if, like me, you self-identify with our community, then happy birthday to you!

TLDR: In 2019 I wrote a post with seven reasons policy professionals should learn to code. Now that AI can write reasonable code, I’ve had a number of people ask me if it’s still worth learning to program. Based on the evidence I’ve seen, AI appears to be a an amplifier, rather than a substitute for human intelligence in the policy space: meaning the people that get the most out of AI are likely to be those that know how to instruct AI to tackle a problem and tell good output from bad. So my answer is yes, learning to code still makes sense. Not necessarily so you can write the code yourself, but so you can make AI a force multiplier for your work.

Note: I’m intending to make updates and refinements to this post over time: both to correct the typos that I haven’t picked up and to incorporate feedback from other coders in the policy community. Feel free to reach out if you have suggestions.

How AI was used: I wrote the first draft of this post. AI was used (minimally) to refine how some ideas were communicated.

Background

In 2019 I wrote a listicle outlining why people working on public policy should learn how to program.

Putting aside that R python was probably the answer to which language to learn, the point of the post was to share some of the practical benefits I’d seen as an applied economist and policy advisor from learning to code.

Despite the listicle being one of my lower effort posts, it seemed to strike a chord: my inbox was quickly flooded with messages thanking me for articulating the benefits others in the public policy community had seen from leveraging code in their work. To my surprise, the article was also used as a reference for several university programming courses, resulting in my inbox still receiving messages from people starting their coding journey.

But a lot has changed since I first wrote the article: the UK voted to leave the European Union; the World Health Organization declared a global pandemic; and 25 movies in the Marvel Cinematic universe were released.

Oh! And a quaint lil’ startup called OpenAI launched something called ChatGPT.

And while there’s certainly a lot that could be said about the unfolding Marvel Universe’s implication for learning to code, this will have to be a subject of a future post. As today, I’d like to attempt to answer a question I’ve been increasingly asked since Artificial Intelligence (AI) tools and/or Large Language Models (LLMs) have made producing code trivial.

It’s sunk costs all the way down

At the outset, I’m clearly not an impartial observer. I was motivated enough to teach myself R, to write the original post and create a dedicated MOOC on using R for policy analysis. I also enjoy coding, so I have a vested interest in answering yes to the question. But, in my defence, these reasons also mean that I’ve thought about it a lot. I’m also a frequent user of AI in my work and find it has greatly expanded the scope and scale of what I can achieve.

This post is therefore my attempt to provide an answer to the question Is it still worth learning to code in 202X? My answer is, yes, but for different reasons to the 2019 listicle.

Doing smart stuff quickly

What got me excited about learning to program when I first started, was its ability to expand the scope and scale of analysis I could do and questions I could answer. In the world of public policy, there’s an almost unlimited number of interesting questions that could be asked, but rarely the time and data to answer them. Sometimes this is simply because the data didn’t exist, which makes learning to code not particularly helpful. However, often there are cases where data does exist, just not in a format that can be readily analysed with the commercial tools made available to you:

• Want to import thousands of individual excel files and merge them into a single dataframe?
• Have a repetitive data cleaning task that would take hours to do using a spreadsheet?
• Want to apply a statistical technique that hasn’t been added to your stats software?

Programming languages makes solving these problems trivial, but with the right prompt, so can AI.

Probability-based slop machines

But let’s be clear: at their core, LLMs are probabilistic slop machines. I don’t mean slop in the pejorative sense, but as a palpable description of what LLMs do: slap together intelligent looking responses based on next-token prediction.[1] I say intelligent looking[2] as LLMs have been trained to produce outputs that look intelligent, rather than to exhibit intelligence in a human sense:

“…AI systems have demonstrated their capability to solve cognitive ability test problems, primarily through guided training (e.g., Zhuo & Kankanhalli, 2020) or programmed approaches to transform problems into algorithmically solvable formats (e.g., Schmidhuber, 2004). While remarkable, it is debatable whether these accomplishments signify intelligence, given that the capabilities of most current AI systems are limited to specific programming and/or training data, without the necessary demonstration of novel problem-solving ability characteristic of human intelligence (Davidson & Downing, 2000; Raaheim & Brun, 1985). Consequently, many AI systems might be more aptly recognised as having the capacity to exhibit artificial achievement or artificial expertise…”

Gignac, G.E. and Szodorai, E.T., 2024. Defining intelligence: Bridging the gap between human and artificial perspectives. Intelligence, 104, p.101832 (link).

But sometimes a little bit of AI slop can be helpful. When producing code, I’ll often get AI to produce the first draft before making tweaks to the approach and style.[3] I use AI to help with writing and research. For this post I asked AI for good search terms to find research exploring the cognitive and economic benefits of learning to code, so I could check my assumptions and think clearly about what has changed since I first wrote the post. I then asked AI for ideas to better communicate key ideas presented in this post.  

“The majority of human production has always been slop. Mediocrity is not a bug of technology; it is the baseline of culture. The canon of art we revere today – those few thousands works in museums and textbooks – is the surviving tip of an immense iceberg of forgotten, derivative, or simply boring creations…”

D’Isa, Francesco (1 December 2025). “The Idea of ‘AI Slop’ Is Slop”. The Philosophical Salon (link)

If you’re wondering why I’m waxing lyrical about intelligence and the mechanics of LLMs, it’s because I think they help decide two things: which problems AI is likely to be good at solving and whether (or when) it’s worth learning to code. After all, if LLMs can already code better and faster than you, why bother learning how to do it yourself? Particularly if the type of intelligence AI outputs can mimic have sufficient crossover with our own to make it practically indistinguishable in the realm of applied economics and public policy.

Widening the scope, scale and quality(?) of analysis

Looking back at my original post, almost all the cited benefits come from programming languages changing the scope, scale and quality of what you can produce as an applied policy analyst. And if I had to put a number on it: of the seven reasons cited in the listicle, six can be achieved without having to write the code yourself.[4] Have some messy data that needs to be cleaned and reshaped before analysis? Share a sample with AI and let it write the code to take care of the data cleaning for you. Have research question you want answered? LLMs can suggest alternative statistical approaches you might not have thought of and write the code to implement them. Need to produce a set of publication ready plots? AI can write the code to do this too in a format that tells the story of your choice.

It would seem to be game over for learning to code and my 2019 post.

Of course, you already know that I don’t believe this but let me tell you why…

The jagged frontier of AI-assisted public policy analysis

Firstly, whether AI is intelligent or merely produces outputs that make it appear to be, early evidence suggests its ability to substitute for human intelligence varies across tasks. Researchers that tested the productivity benefits of AI termed this as a jagged frontier, that the users had to carefully navigate to benefit from using AI in their work.[5] It’s also not clear whether improving model performance is likely to expand the frontier equally across tasks, with many of the benchmark scores not necessarily generalizing outside test questions and being fragile to irrelevant and/or incorrect context.[6] Finally, at the time of writing this post, the types of tasks that frontier AI models had the hardest time completing successfully, looked strikingly similar to the problems encountered in public policy: messy, interdependent and heavily reliant on context.[7]

Secondly, while there’s solid evidence suggesting that AI can lead to productivity gains for knowledge workers, this isn’t automatic and appears to lean on the taste, judgement and experience of the user. [8] This was demonstrated by a field experiment in Kenya,[9] where a little over 600 entrepreneurs were randomly provided access to an AI-powered business assistant. The performance of those with the assistant was then tracked over time and compared to those without the AI assistant. The result: top performers gained around 15 percent from having access to AI, while weaker ones lost around 10. Entrepreneurs that benefited the most from AI appeared to be better at identifying and acting on contextually appropriate advice, while poor performers tended to adopt generic strategies proposed by the tool.

A much smaller experimental study closer to the world of coding pointed to a similar idea,  with both the 2025 (n=16) and 2026 (n=57) replication suggesting the productivity returns to using AI are ambiguous and vary based on the experience of the user.[10] Both versions of the study took the same basic form: experienced developers were randomly assigned to work with AI or not on their own repositories and the performance benefits were compared between groups. In the 2025 study, the authors found the performance benefits were on average negative: experienced developers took longer to complete the same tasks and overestimated the time saved.  In the 2026 replication,[11] the authors found performance increased on average, but not for all developers.[12] Although the authors note their estimates are likely biased downwards due to difficulties in recruiting developers willing not to use AI, they paint the same basic picture: the benefits of using AI are not unambiguously positive and depend on the skills and experience of the user.[13]  

“…The right analogy for AI is not humans, but an alien intelligence with a distinct set of capabilities and limitations. Just because it exceeds human ability at one task doesn’t mean it can do all related work at human level. Although AIs and humans can perform some similar tasks, the underlying “cognitive” processes are fundamentally different…”

Mollick, E. (2024, May 12). Superhuman? What does it mean for AI to be better than a human? And how can we tell? One Useful Thing (link).

There’s no accounting for taste

To be clear, my point here isn’t to suggest AI tools aren’t helpful, just that the size of the benefit depends on the judgement, expertise and taste of the user. AI certainly can suggest a variety of statistical techniques for answering your research question, but are they the right ones?[14] Sure, AI can quickly produce code for twenty reasonable looking plots, but do you understand the precise changes needed to accurately present your results and apply your corporate style guidelines and colours? Certainly, AI might simplify data cleaning, but does the approach make the right assumptions given the context, your research question and how the data was collected? And, while it might be true that research built on AI-generated code is more reproducible than a spreadsheet with a tangled collection of formulas, it may just be encoding hidden errors, the wrong methodology and creating a new source of technical debt you and your team will have to pay in the future.

In my original post I was excited at programming expanding the scale, scope and quality of analysis I could do. AI can help with the first two, but its influence on the latter seems to depend on the user. AI may be able to slap together a solid draft that gets you 90 percent of the way, which for some tasks is all you need, but deciding what good looks like and instructing AI how to help relies on building sufficient judgement, taste and expertise. 

Vibing fast and slow

The need to develop judgement, expertise and professional taste is also not unique to coding. But like musicians with sheet music, a mathematician with an equation and an economist with a set of national account statistics, anyone working with code is likely to perform better with some understanding of how it works, even with the help of AI. In the words of Andrew Ng:

“…people who understand the language of software through their knowledge of coding can tell an LLM or an AI-enabled IDE what they want much more precisely, and get much better results.”

And once again, there’s research that backs this up. For instance, pre-AI a meta-analysis of the effects of learning to code points to benefits to creativity and metacognition i.e. the processes underlying the monitoring, adaptation, evaluation, and planning of thinking and behaviour,[15] both which are plausibly linked to the ability to properly prompt and manage AI to perform tasks. A small observational study[16]  (n=100) seems to support this idea, with their analysis finding a positive association between vibe-coding efficacy and computer science (CS) domain expertise. Although the exact questions used to evaluate CS expertise aren’t explicitly outlined, their cited source points to them covering the same thinking and skills that learning to program helps to build,[17] with the authors concluding that effectively vibe-coding likely draws on problem decomposition and algorithmic thinking skills.

The book Messy Jobs also makes this point by noting that by cheapening the cost of cognition, AI simultaneously raises the value of judgement and taste.[18] The consequence of this for vibe-coding is that both instructing AI how to solve a problem and being able to judge the quality of outputs relies on knowing how code works, otherwise you might just be doing stupid things faster, in the words of Hadley Wickham:

“…LLMs and AI are, I think, a massive amplifier. You can certainly use them to amplify the worst aspects of yourself, but you can also use them to amplify the best aspects…”

Wickham, H. (2026, July 25). y code when ai? A talk on AI and coding. Tidy Design. https://tidydesign.substack.com/p/y-code-when-ai

Summing up

In my 2019 listicle, I argued that learning to code was worth it so you can answer a wider variety of questions more quickly and at a higher standard to traditional point-and-click software. Seven years later and AI tools can produce code for a wide variety of tasks more quickly than an experienced developer. But that doesn’t make what it produces good.

None of this means that I think you should avoid using AI before having mastered how to code, just that to be a force-multiplier it helps to know what good code looks like. Treat it as a brainy intern: motivated, well read, and in need of supervision. Give it clear instructions and it might come back with something useful. Don’t and its code might set the office microwave on fire. Learning to code can help you tell the difference between the two.

That matters because you’re accountable for the result as the human given the the instructions. Use AI to decide how to prioritize which grants to cancel and you might cause unnecessary harm,[19] trust its recommendations without validating them and you might just erroneously reinforce cultural stereotypes[20], and have it author a report for you, and it might just cite hallucinated sources.[21] If used well, AI can help you reach your destination quicker and expand your professional frontier. If used poorly, it might just having you producing more errors at scale.

I don’t know how long this holds. The studies here are small and not focused directly on the tasks we do in public policy. My best guess is that the jagged frontier will continue to expand, but unevenly, which will continue to place a premium on professional expertise and taste for tasks AI continues to struggle with.

So: should you still learn to code in 202X? Yes, but not for the reason I gave in 2019. Learn to code so you can tell good code from bad and steer the tool in the right direction.

[1]The push-back I’ve seen to this description of how LLMs work doesn’t change the core point. LLMs are trained to probabilistically predict the next token to produce intelligent-looking outputs at scale.

[2] LLMs are trained on the finished product of human intelligence, such as published research, books and social media posts, but not the thinking that went into producing them. Newer models that ‘think’ before responding to a prompt are meant as a partial workaround to this, but don’t change the basic point that LLMs are trained to mimic articulated intelligence, rather than the process of producing. 

[3] A practice that might have the opposite result if you’re not careful, see: Goldsmith-Pinkham, P. (2026, September 10). How do we know if writing is AI? Some initial analyses of AI’s impact on economics writing. A Causal Affair (link) 

[4] The one that may no longer be true is the trend towards open-source statistics. Certainly, I’ve been experimenting with ways to use AI to better extract government data (see example), but it’s entirely possible LLMs will make the internet a more difficult resource to navigate as text can be generated at larger scales.

[5] Dell’Acqua, F., McFowland III, E., Mollick, E., Lifshitz, H., Kellogg, K.C., Rajendran, S., Krayer, L., Candelon, F. and Lakhani, K.R., 2026. Navigating the jagged technological frontier: Field experimental evidence of the effects of artificial intelligence on knowledge worker productivity and quality. Organization Science, 37(2), pp.403-423.

[6] Akimitsu, P., 2026. Wrong and More Confident: A Field Experiment on Large Language Models Taking a Graduate Economics Exam. arXiv preprint arXiv:2607.23424.

[7] METR. (2026, May). Task-completion time horizons of frontier AI models. https://metr.org/time-horizons/

[8] Ethan Mollick uses the term Centaurs and Cyborgs todescribe the two most common approaches used by  workers that were able to navigate this frontier. In both cases, using AI effectively requires intelligently and selectively integrating the tools into their work without ‘falling asleep at the wheel’ (link).

[9] Otis, N.G., Clarke, R., Delecourt, S., Holtz, D. and Koning, R., 2023. The Uneven Impact of Generative AI on Entrepreneurial Performance: Evidence from a Field Experiment in Kenya (No. 24-042). Harvard Business School Working Paper.

[10] Becker, J., Rush, N., Barnes, B., & Rein, D. (2025, July 10). Measuring the impact of early-2025 AI on experienced open-source developer productivity. METR. https://metr.org/blog/2025-07-10-early-2025-ai-experienced-os-dev-study/

[11] Becker, J., Rush, N., Cunningham, T., Rein, D., & Mahamud, K. (2026, February 24). We are changing our developer productivity experiment design. METR. https://metr.org/blog/2026-02-24-uplift-update/

[12] The authors note it had become increasingly in their 2026 follow-up study to recruit developers willing to work without AI, which likely resulted in estimates of AI-assisted speedup being biased downward.

[13] Although the authors recommend their results are interpreted with caution, the 2025 study not only found that developers using AI were slower, but that more experienced developers tended to overestimate the time they would save from using AI.

[14]Like humans, LLMs exhibit systematic biases in the way they work and evaluate problems. This might come in subtle forms that influence outputs subtly, such as assuming doctors are more likely to be men than women. Or explicit, such as producing outputs that reflect how prominent a solution is in its training set, rather than necessarily being appropriate for the task at hand (link).  

[15] Scherer, R., Siddiq, F. and Sánchez Viveros, B., 2019. The cognitive benefits of learning computer programming: A meta-analysis of transfer effects. Journal of Educational Psychology, 111(5), p.764.

[16] Thorgeirsson, S., Weidmann, T.B. and Su, Z., 2026, April. Computer science achievement and writing skills predict vibe coding proficiency. In Proceedings of the 2026 CHI Conference on Human Factors in Computing Systems (pp. 1-17).

[17] Specifically, code tracing, code completion, conditionals, loops and logical operators, see: Parker, M.C., Solomon, A., Pritchett, B., Illingworth, D.A., Marguilieux, L.E. and Guzdial, M., 2018, August. Socioeconomic status and computer science achievement: Spatial ability as a mediating variable in a novel model of understanding. In Proceedings of the 2018 ACM Conference on International Computing Education Research (pp. 97-105).

[18] Garicano, L, Li, J & Wu, Y 2026, Messy Jobs: The Work That AI Cannot Reach, Upriver Press

[19] Alibašić, H., 2026. Introduction: The Imperative of Hybrid Intelligence in Digital Governance. In Hybrid Intelligence for Effective Digital Governance: AI in Administration (pp. 3-48). Cham: Springer Nature Switzerland.

[20]Saar Alon-Barkat, Madalina Busuioc, Human–AI Interactions in Public Sector Decision Making: “Automation Bias” and “Selective Adherence” to Algorithmic Advice, Journal of Public Administration Research and Theory, Volume 33, Issue 1, January 2023, Pages 153–169, https://doi.org/10.1093/jopart/muac007 

[21] Karp, P. (2025, November 6). AI-tainted Deloitte report was worse than previously thought. Australian Financial Review. https://www.afr.com/politics/ai-tainted-deloitte-report-was-worse-than-previously-thought-20251106-p5n863

The post X Reasons for policy professionals to get into programming in 202X appeared first on Giles.

Dear rOpenSci friends, it’s time for our monthly news roundup! You can read this post on our blog. Now let’s dive into the activity at and around rOpenSci!

rOpenSci HQ

Yani and Mark at “Open Communities in the age of AI”

rOpenSci community manager Yani and software-review lead Mark will participate in the upcoming Openscapes community call on October 1st (this Thursday!) at 9:30AM PT (16:30 UTC) entitled “Open Communities in the Age of AI”, together with Mara Averick, Senior Developer Advocate at Quansight and Hadley Wickham, Chief Scientist at Posit.

Event page, including link for free registration.

Champions Program update

Our current 2026–2027 cohort has completed the training phase and is now focused on developing their individual projects with the support of their mentors, as well as on their outreach activities.

Meanwhile, the 2025–2026 cohort will wrap up their journey with a closing community call, Más Allá del Código: muestra abierta de los proyectos de nuestros campeon(a|e)s (Beyond the Code: an open showcase of our Champions’ Projects). At the call, the Champions will present their final projects, from creating, improving, and reviewing R packages to the outreach work that brought those projects to their communities. They’ll also share what their projects led to, including scholarships, new collaborations, and professional growth. The event is in Spanish. The call will take place on Monday, October 12 at 15:00 UTC. Come get inspired, ask your questions live, and help us celebrate our growing open-source community!

We’re still celebrating our 15th anniversary!

In July, we started to share stories from members of our community about their experiences with rOpenSci. Our second story features Yi-Chin Sunny Tseng and her connection with rOpenSci. Read it on our blog: Happy Birthday rOpenSci — My Journey from First-Time Developer to Current Opportunities Stay tuned for more stories from our community as we continue celebrating 15 years of rOpenSci!

rOpenSci at LatinR 2026 in Medellín, Colombia. See you there!

Registration for LatinR 2026 (November 11–13, Universidad de Antioquia, Medellín) is now open, and rOpenSci will have a strong presence.

Jeroen will give a talk on R-Universe. Yani will lead a workshop on R package development and give a talk on the Champions Program’s open curriculum. Nic Crane is one of the conference’s keynote speakers. Evelia Lorena Coss Navarrete, one of our Champions, will talk about the package she developed during the program. Several other rOpenSci folks — Natalia Da Silva, Luis Verde, and Francisco Cardozo, among others — will also be there. Jeroen, Francisco, and Yani will take part in the hackathon during the conference.

Join us in Medellín!

Coworking

Read all about coworking!

• Tuesday October 6th, 09:00 Americas Pacific (16:00 UTC) “Writing Tests & Testing in R”, with Yanina Bellini Saibene and co-host Olivier Leroy.
  • Explore how to write tests for R and add some tests to your work or packages
  • Meet co-host, Olivier Leroy, and chat about testing
• Tuesday November 3rd, 09:00 Australia Western (01:00 UTC) “Climate Science in R”, with Steffi LaZerte and co-host Elio Campitelli.
  • Explore how R is used to study the climate
  • Meet co-host, Elio Campitelli, and discuss Climate Science in R
• Tuesday December 8th^*, 14:00 Europe Central (12:00 UTC) “Code Linting in R”, with Steffi LaZerte and co-host Etienne Bacher.
  • Read up on Code Linting and apply some linters to your R code
  • Meet co-host, Etienne Bacher, and discuss code linting in general, or flir and Jarl in particular
    * Note that December coworking is a week later than usual

And remember, you can always cowork independently on work related to R, work on packages that tend to be neglected, or work on what ever you need to get done!

Software

The following two packages recently became a part of our software suite:

• ciecl, developed by Rodolfo Tasso Suazo: Tools for working with the International Classification of Diseases (ICD-10 Chile official MINSAL/DEIS v2018). Includes optimized SQL search with SQLite, fuzzy matching of medical terms (Jaro-Winkler), Charlson and Elixhauser comorbidity calculation, WHO ICD-11 API integration, and hierarchical code validation. Data from Centro FIC Chile DEIS https://deis.minsal.cl/centrofic/. It has been reviewed by Maëlle Salmon and Yanina Bellini.

• brapiR2, developed by Joash Joshua Ayo: Provides pipe-friendly, stateless read access to the Breeding API (BrAPI) v2.1 specification, an open community standard for plant breeding data interchange maintained by the BrAPI project https://brapi.org. Wraps 32 of the 37 BrAPI v2.1 entities across all four modules, Core, Germplasm, Phenotyping, and Genotyping, covering 49 of the specifications 138 retrieval (GET and search) endpoints and returning tidy tibbles ready for analysis. Write and update endpoints are out of scope by design. Features include automatic pagination, async search handling, response caching, parallel batch fetching, and convenience functions for genomic selection workflows (e.g. dosage matrix extraction). Designed for plant breeders and bioinformaticians who need programmatic access to plant breeding databases that implement the BrAPI’ v2 specification. It has been reviewed by David Waring and Jenna Hershberger.

Discover more packages, read more about Software Peer Review.

New versions

The following twenty-two packages have had an update since the last newsletter: visdat (v0.6.1), nycOpenData (v0.2.3), ruODK (v1.6.0), stats19 (v4.1.0), npi (v0.3.1), promoutils (v0.7.0), brapiR2 (v0.2.0), ciecl (v1.0.0), bibtex (v0.5.3), ckanr (v0.9.0), osmapiR (v0.2.6), distionary (v0.2.0), GLMMcosinor (v0.2.2), reviser (v0.3.1), autotest (v0.2), occCite (v0.6.3), readODS (v2.3.6), osmdata (v0.4.1), goodpractice (v1.2.0), galamm (v0.4.1), RAMEN (v2.1.2), and EDIutils (v3.0.1).

Software Peer Review

There are seventeen recently closed and active submissions and 4 submissions on hold. Issues are at different stages:

• Two at ‘6/approved’:

  • brapiR2, A Tidyverse-Native Client for the BrAPI v2 (Breeding API) Specification. Submitted by Joash Joshua Ayo.

  • ciecl, International Classification of Diseases ICD-10/ICD-11 for Chile. Submitted by Rodolfo Tasso.

• Three at ‘5/awaiting-reviewer(s)-response’:

  • ibger, Access the IBGE Aggregate Data API from R. Submitted by Andre Leite Wanderley.

  • RAQSAPI, A Simple Interface to the US EPA Air Quality System Data Mart API. Submitted by mccroweyclinton-EPA.

  • coevolve, Fit Bayesian Generalized Dynamic Phylogenetic Models using Stan. Submitted by Scott Claessens. (Stats).

• One at ‘4/review(s)-in-awaiting-changes’:

  • rcrisp, Automate the Delineation of Urban River Spaces. Submitted by Claudiu Forgaci. (Stats).

• Four at ‘3/reviewer(s)-assigned’:

  • camtrapReport, Camera-Trap Report Generator. Submitted by Elham Ebrahimi.

  • nert, Curated Access to TERN Environmental Raster Data. Submitted by Max Moldovan.

  • rfastlowess, High-Performance LOWESS Smoothing for R. Submitted by Amir Valizadeh. (Stats).

  • camtrapReport, Camera-Trap Report Generator. Submitted by Elham Ebrahimi.

• Two at ‘2/seeking-reviewer(s)’:

  • grumpy, Read NumPy .npy and .npz Files. Submitted by Hugo Gruson.

  • tezr, Access Thesis Metadata from Turkiye’s National Thesis Center. Submitted by Emrah Er.

• Five at ‘1/editor-checks’:

  • ocean3d, Three-Dimensional Marine Spatial Analysis. Submitted by Jay Matsushiba.

  • corila, Sparse modelling with grouped and correlated features allowing for privileged information. Submitted by Armin Rauschenberger. (Stats).

  • tarpolyglot, Run Python, Julia, Rust, and C++ Inside targets Pipeline Steps. Submitted by Pierre9344.

  • OptSurvCutR, Optimal Survival Cut-Point Discovery for Time-to-Event Analysis with OptSurvCutR. Submitted by Payton Yau. (Stats).

  • HydraR, Stateful Agentic Orchestration for Scientific Reproducibility. Submitted by Ignatius Pang.

Find out more about Software Peer Review and how to get involved.

On the blog

• Happy Birthday rOpenSci — My Journey from First-Time Developer to Current Opportunities by Yi-Chin Sunny Tseng. Yi-Chin Sunny Tseng shares how the rOpenSci Champions Program helped her grow from a first-time R package developer into an open science contributor, creating tools for biodiversity research and discovering new opportunities along the way.

Calls for contributions

Calls for maintainers

If you’re interested in maintaining any of the R packages below, you might enjoy reading our blog post What Does It Mean to Maintain a Package?.

• charlatan, create fake data in R. Issue for volunteering.

Calls for contributions

Refer to our help wanted page – before opening a PR, we recommend asking in the issue whether help is still needed.

Package development corner

Some useful information for R package developers.

{covr2gh}: Coverage summary as GitHub comment

The covr2gh package by Dragoș Moldovan-Grünfeld provides an automated way to summarise the impact of a pull request (PR) on test coverage directly within GitHub. See it in action in a PR in the tfrmt R package.

Deprecation messages for package data

Hugo Gruson wrote an exhaustive post about the deprecation of data in an R package. The post features the delayedAssign() function.

A refactoring story featuring people

Athanasia Mo Mowinckel published “Why ggseg Atlases Became Function Calls”, where she explains how she made data into objects into a package in order to allow re-exporting them. She furthermore tells how she got to that conclusion by trying out different solutions and discussing with other package developers in the rOpenSci Slack workspace.

usethis 3.2.2

The usethis package was updated on CRAN. The changelog features many quality of life improvements around Git and GitHub, and a new use_readme_qmd() function for drafting a README in the Quarto format.

Last words

Thanks for reading! If you want to get involved with rOpenSci, check out our Contributing Guide. This guide will help direct you to the right place, whether you want to make code contributions, non-code contributions, or contribute in other ways such as through sharing use cases. You can also support our work through donations.

If you haven’t subscribed to our newsletter yet, you can do so though our signup form. Until it’s time for our next newsletter, you can keep in touch with us through our website, Mastodon, or LinkedIn. See you soon!

This is a guest post by Matt Denwood.

After a longer-than-anticipated delay, the macOS installer for JAGS 5.0.0 is now available via SourceForge. There are some big changes for how the official binaries of JAGS 5.x are installed on macOS, which is the topic covered by this news post.

Multi-threading and BLAS linkage

JAGS now has multi-threading support to allow parallel execution of independent chains, which uses Apple’s Grand Central Dispatch for the macOS binaries. This is the default build, as most users will want to benefit from the parallel computation. However, a single-threaded build is also provided as part of the installer, which has the same behaviour as JAGS 4.3.2. A third build is also provided within the installer, which uses single threading along with Netlib’s reference implementation of BLAS and LAPACK rather than Apple’s vecLib (Accelerate framework). By default, all three of these builds are installed so that the user can later switch between them;  however, the installer can be customised to deselect builds that are not needed.  For any users wishing to avoid the larger download and installation footprint (e.g. build/test machines), we have also provided single-build installers via SourceForge. For discussion on the difference between BLAS implementations, see the relevant section of the R for macOS FAQ. It is also possible to build JAGS using OpenMP on macOS – this is not provided as part of the official macOS installer, but the necessary compilation instructions will be added to the macOS section of the JAGS installation manual in due course.

Installation location

1. Although there is no standardised location for installing system-level software on macOS, the /opt directory is now more frequently used (by e.g. homebrew and some elements of R) than /usr/local. Part of the reason for this is that installation under /opt/jags ensures a clean separation between software relating to JAGS and other software that you may have installed under other subdirectories of /opt.
2. The installer now contains multiple builds of JAGS as well as utilities that are specific to macOS (see below) – the versioned directory structure allows the different builds and utilities to coexist.
3. Having a versioned directory structure also allows different major versions of JAGS to be installed on the same system – we hope that this will ease the transition from JAGS 4.3.2 to JAGS 5.0.0 on CRAN’s build systems.

The installer also contains a build of JAGS 4.3.2 that installs under the same /opt/jags directory, with symbolic links to the older /usr/local directory that allow the CRAN build of rjags-4 to continue to function as expected. These symlinks also allow JAGS 5.x to be accessed from the Terminal (and by pkg-config) without updating the PATH variable, although we recommend adding the new directory to the PATH explicitly using (for example): echo 'export PATH="/opt/jags/bin:$PATH"' >> ~/.zshrc

JAGS utilities

As well as JAGS itself, the macOS installer provides the following three utilities:

• jags-version – this allows the active build of JAGS to be switched, both for the command-line version and for rjags.  Switching between major versions of JAGS requires re-installation of the rjags package with matching version, as well as re-compilation of the runjags package. Switching between builds of JAGS with the same major version does not require re-installation of any R package (but may require R to be restarted, if the required packages are loaded).
• jags-uninstall – this allows for complete removal of any JAGS installations detected either under the new /opt/jags directory or the legacy /usr/local directory.
• pkgconf-lite – as macOS does not include pkg-config by default, this small utility (provided by pkgconf) is installed in order to facilitate compilation of rjags on vanilla macOS installations.

See the provided man pages for jags-version and jags-uninstall for more information on these utilities.  Note that these utilities are not installed if using the single-build installers.

Getting help

If you encounter any issues with the macOS installer or utilities, please either post to the macOS-specific JAGS installation thread on the JAGS discussion forum. If you have comments or suggestions regarding the macOS utilities or build process for JAGS then please feel free to open an issue on the GitHub repo for the macOS build scripts.

I had a lot of fun this morning, brainstorming and assembling this Analytics Pipeline for Dashboards with Claude (and yes it takes much, much more than only 10 prompts).

The philosophy is borrowed from Observable Framework’s data loader, but implemented from scratch with voluntarily opinionated choices of libraries for Python, R and Javascript (and no Markdown): pre-compute data with Python + Polars (Polars is my new crush) or R + dplyr, publish it as static files, and explore them with JavaScript in the browser (Arquero for data wrangling, Observable Plot or Highcharts for charts)

The repository is available on GitHub and here’s the quick start guide to run it locally (you need to have Python and R installed on your machine):

make                        # list all targets
uv venv venv                # or: make venv (python -m venv venv without uv)
source venv/bin/activate    # Windows: venv\Scripts\activate
make install                # Polars into venv/, dplyr into .uvr/library/
make dev                    # starts http://localhost:3000 and opens your browser

The repository’s README file is also informative and contains a few more details about the pipeline.

There are ten public training courses left in 2026, running through October and November. If you’ve been meaning to learn Python, build your first Shiny app, or get to grips with machine learning in R, there’s still time to do it before the end of the year.

All courses run online and are delivered live by data scientists and engineers who work on client projects. Every participant gets course notes and scripts, live demonstrations, and hands-on exercises throughout.

What’s On This Autumn

Python in October

October is a good month to pick up Python. Introduction to Python runs on 12th–13th October, Programming with Python follows a week later, and Data Visualisation with Python runs straight after that. You can take all three in turn to go from your first lines of Python to building your own functions and charts in about two weeks.

Shiny

Introduction to Shiny runs on 5th–6th October and takes you from an empty script to a working interactive web app in R. If you already build Shiny apps, Advanced Concepts in Shiny on 23rd–24th November covers how to write maintainable code, build robust apps that handle unexpected user input, and get more out of reactive programming.

Visualisation and Statistics with R

Data Visualisation with ggplot2 runs on 2nd–3rd November. It’s a good fit for anyone who uses R already and wants to produce clear, publication-quality charts. Statistical Modelling with R runs on 16th–17th November and covers hypothesis testing, regression, clustering and principal components analysis.

Machine Learning with Tidymodels

Machine Learning with Tidymodels runs on 4th–5th November, and Advanced Machine Learning with Tidymodels follows on 23rd–24th November. Taking both gives you a complete machine learning workflow in R, from pre-processing and model fitting through to tuning.

Git

Introduction to Git runs on the mornings of 16th–17th November. It covers why version control matters, how to use it properly, and how to start collaborating with others using Git.

Git and Advanced Machine Learning with Tidymodels run in the mornings, so you can take them alongside an afternoon course in the same week.

Autumn 2026 Schedule

Each course runs over two sessions of three and a half hours. Afternoon sessions run from 13:30 to 17:00 and morning sessions from 09:00 to 12:30 (UK time). Click a course to see the full outline and book.

Booking closes one week before each course starts, so Introduction to Shiny closes on 28th September.

Early bird pricing is still available on some November courses: until 4th October for Introduction to Git and Statistical Modelling with R, and until 11th October for both Advanced courses.

Planning Ahead? Book Now for 2027

If the autumn is too busy, you can already book for next year. Our public schedule for January to June 2027 is live, and every course on it is currently at the early bird price.

It follows the same pattern as this year. Introduction to R and Data Wrangling in the Tidyverse run in January and again in April, and the Python courses run in February and May. Machine learning with Tidymodels runs in March.

It also includes Introduction to Bayesian Inference using RStan in January. It runs over four afternoons and covers writing Stan programs for a range of statistical models and checking MCMC diagnostics in R.

Early bird pricing ends about six weeks before each course, so for the January courses that means booking by late November or early December.

We’ll add dates for the second half of 2027 later in the year.

Why Train with Jumping Rivers

Our trainers are practising data scientists and engineers. The examples and exercises in every course come from the work we do with clients, so you learn how these tools are used in real projects.

Every participant receives:

• Comprehensive PDF notes and scripts to keep after the course
• Live demonstrations and hands-on exercises throughout
• Direct access to a trainer for questions

We have delivered over 1,000 courses to organisations including NHS Scotland, Shell, Wessex Water and the Royal Statistical Society.

Additional Perks

We run free webinars throughout the year. Attending a Jumping Rivers webinar gives you:

• Early exposure to new topics in data science and analytics
• Up to 20% off training courses
• Up to 20% off Jumping Rivers conferences

Register for upcoming webinars here.

Training for Teams

We also run in-house training for organisations that want to develop their teams. We tailor each course to your workflows, tools and experience levels, and we can deliver any course in our full catalogue, not just the ones on the public schedule. Group bookings and returning clients get a discount.

To discuss options, email training@jumpingrivers.com.

Book Your Place

View the full schedule and book at jumpingrivers.com/training/public. Places on public courses are limited, so book early if you have a particular date in mind.

For updates and revisions to this article, see the original post

Before the main content:

• I am creating an R Community on Google Groups. You can join the group using this form.
• With all the funding cuts to PhD studies, I’d appreciate if you can donate here.

In Economics we often use reduced-form models, which are simplified ways to describe a process. I tend to see the criticism to this as unfair, as in physics they do the same, and so many models in physics assume “spherical objects in a vacuum” as a simplification to make real-world problems mathematically solvable.

With reduced form models we are not interested in how faithful those ressemble reality but how useful these are. As an example in support of reduced form models, consider the two extreme cases of models to answer “how to get from 37 Old Queen St to St. James Park underground station?”:

1. Person 1 draws a quick sketch like this on an old receipt

      _ Embassy
     |   
  ---
 |
Underground

1. Person 2 built a very detailed 1:10000 architectural model of St. James neighbourhood and shows you the tiny streets, cars, etc.

For the question, person 1 sketch is sufficient.

The Armington model of trade sits between these two extremes. It considers a world of \(N\) countries, each country \(i\) has a fixed labour endowment \(w_i\) and there is an endogenous labour wage. Income is given by \(Y_i = w_i L_i\). Trade flows are given by \(X_{ij}\) (and so the sales share is \(\gamma_{ij} = X_{ij} / Y_i\)). There is a productivity shift \(\x_i\), a trade cost $_{ij}, and an elasticity of substitution \(\varepsilon\).

The strong assumption in this model is a known and homogenous elasticity of substitution. In this model countries can substitute between Apple and Android phones even when many would argue that those phones compete in different markets.

The model imposes a market clearing condition \(w_i L_i = \sum_{j=1}^N \lambda_{ij} w_j L_j\).

There are many values \(\lambda_{12}, \lambda_{13}, \ldots, \lambda_{1N}\) that balance the market clearing condition but in this model the solution is a gravity-type equation that leads to ${ij} = (i / ({ij} w_i)^{}) / sum{l=1}^N (l / ({lj} w_l)^{}). More about gravity-type equations can be read in Gravity for Undergrads by Prof. Dr. Yoto V. Yotov.

Consider, for example, the EU increase on low value imports from the UK described [https://hboltd.co.uk/uk-eu-shipping-changes-2026-customs-overhaul/]. In this model, instead of a product-by-product specific tariff have an initial UK-EU average trade cost of \(\tau_{ij}\) that increased to \(\tau_{ij}’\). The hat form of this is \(\hat{\tau_{ij}} = \tau_{ij}’ / \tau_{ij}\).

After the trade cost change the market clearing would be \(w_i’ L_i’ = \sum_{j=1}^N \lambda_{ij}’ w_j’ L_j’\), and then

\[ \hat{w_i} \hat{L_i} = \frac{\sum_{j=1}^N \lambda_{ij}’ w_j’ L_j’}{\sum_{j=1}^N \lambda_{ij} w_j L_j} \]

\[ \hat{w_i} \hat{L_i} = \frac{\sum_{j=1}^N X_{ij}’}{w_i L_i} \]

\[ \hat{w_i} \hat{L_i} = \sum_{j=1}^N \gamma_{ij} \hat{X_{ij}} \]

Similarly, $ = .

The new hat equations are useful, for example, to determine the effect of the change to the trade cost over wages \(\hat{w_i}\). The hat market clearing means that \(\hat{w_i} = \sum_{j=1}^N \gamma_{ij} \hat{X_{ij}}\) (labour is exogenous), and therefore \(\hat{w_i} = \sum_{j=1}^N \gamma_{ij} \hat{\lambda_{ij}} \hat{w_j}\).

Merging the hat market clearing with the hat gravity equation, we get

\[ \hat{w_i} = \sum_{j=1}^N \frac{\gamma_{ij} \hat{w_j} / (\hat{\tau_{ij}} \hat{w_i})^{\varepsilon}}{\sum_{l=1}^N \gamma_{lj} \hat{w_l} / (\hat{\tau_{lj}} \hat{w_l})^{\varepsilon}} \]

From the initial equilibrium shares \(\lambda_{ij}\) and \(\gamma_{ij}\) and a fixed \(\varepsilon\) (for example, estimated using USITC data), we have a system of \(N \times N\) equations that we can solve to get \(\hat{w_j}\) and then recover \(\hat{\lambda_{ij}}\). In other words, we simplified an identification problem by solving the differences by resorting on a sufficient statistic instead of estimating all trade costs.

Executive Summary

MicroStrategy (MSTR) has fundamentally transitioned from a traditional enterprise software firm into an equity-based Bitcoin (BTC) holding vehicle and treasury operation. Evaluating MSTR using conventional corporate finance metrics (such as Price-to-Earnings or EBITDA multiples) fails to capture the core driver of its equity valuation: the dynamic Net Asset Value (NAV) premium driven by programmatic Bitcoin treasury expansion.

This technical article explores an end-to-end, reproducible quantitative pipeline built in R. Inspired by recent advancements in automated SEC parsing tools—specifically the secfile architecture highlighted in Interactive Brokers Quant Blog—this pipeline ingests real-time SEC EDGAR filings, dynamically extracts digital asset balance sheet facts, aligns mixed-frequency financial and market data, and models MSTR daily equity log returns using an out-of-sample 82.08% R-Squared structural tidymodels framework.

1. The SEC EDGAR Challenge & The secfile Paradigm

Unstructured Financial Data vs. Machine-Readable XBRL

Historically, extraction of balance sheet metrics directly from SEC filings (Forms 10-K and 10-Q) required brittle HTML regex scraping or costly third-party commercial APIs. Corporate disclosures often vary across reporting periods, particularly when handling emerging asset classes like digital currencies. MicroStrategy’s XBRL taxonomy has evolved, tagging Bitcoin holdings under varying terms such as DigitalAssets, CryptocurrencyHoldings, or within broader Assets line items.

As demonstrated in the Interactive Brokers Quant Blog overview of secfile, direct ingestion of SEC EDGAR financial facts via compliant XBRL parsing solves three primary institutional hurdles:

1. Regulatory Compliance & Transparency: Enforcing strict HTTP User-Agent headers compliant with SEC EDGAR access policies ensures uninterrupted data pipeline ingestion.
2. Dynamic Taxonomy Mapping: By searching the taxonomy dictionary programmatically, the pipeline automatically handles changes in reporting terminology without hardcoded column dependencies.
3. Point-in-Time Alignment: Raw SEC filings record historical submission dates, preventing look-ahead bias when backtesting event-driven trading strategies against historical market prices.

Using tidy evaluation via the rlang package, the pipeline dynamically isolates balance sheet events safely without breaking execution when SEC XBRL tags undergo structural revisions.

2. Mathematical Structure & Structural Economics

The core hypothesis of this quantitative framework is that MSTR daily return dynamics are governed by two orthogonal components:

1. Systematic Asset Return Co-movement: Direct market return exposure to spot Bitcoin returns.
2. Balance Sheet Shock Factors: Discontinuous structural changes in treasury balance sheet holdings resulting from capital raises or debt-funded BTC acquisitions.

Daily Log Returns

Continuous log returns for asset prices are calculated by taking the natural logarithm of the current price divided by the previous day’s price. Taking natural logarithms guarantees additivity over time horizons and prevents non-negative boundary issues inherent to simple percentage returns.

Balance Sheet Shock Factor

To quantify treasury expansion independently of market price fluctuation, the Balance Sheet Shock metric measures the logarithmic growth rate of total Bitcoins held on MicroStrategy’s balance sheet.

Because corporate filings occur at discrete quarterly intervals, daily balance sheet values are carried forward using Last Observation Carried Forward via the zoo package. Consequently, the balance sheet shock value equals zero on non-reporting days, triggering discrete structural impulses only on filings or event reporting dates.

Econometric Specification

The structural linear model specification expresses the daily log return of MSTR as a function of three main components:

• The intercept term, representing base drift.
• The spot Bitcoin log return, weighted by the Equity-to-BTC Beta Elasticity coefficient (measuring leverage and NAV premium response).
• The Balance Sheet Shock factor, weighted by the Treasury Expansion Impact coefficient.
• A residual error term capturing unexplained market noise.

3. Pipeline Architecture & Package Integration

The pipeline follows a modern, production-grade functional programming architecture using the tidyverse ecosystem:

• secfile: High-throughput scraper and XBRL parser for SEC EDGAR datasets.
• tidyquant: Financial wrapper connecting quantmod and PerformanceAnalytics packages into tidy data frames.
• timetk: Provides time-series splitting for chronological out-of-sample data partitioning without temporal leakage.
• tidymodels: Modular modeling framework combining recipe creation, parsnip engine specification, and workflow execution.
• plotly & ggtext: Interactive visualization engine supporting HTML/CSS styled markdown tooltips and dynamic price projection bands.
• zoo: Handles non-homogeneous time-series alignment via Last Observation Carried Forward (na.locf).

4. Chronological Splitting & Model Training

Financial time series violate the Independent and Identically Distributed (i.i.d.) assumption of standard k-fold cross-validation due to autocorrelation. To preserve time order, we employ strict time-based windowing via timetk.

Using recipes and workflows, we define the feature roles and couple them with a native Ordinary Least Squares (OLS) estimation engine. This guarantees clean execution without data leakage across the training split and testing horizon.

5. Empirical Results & Performance Evaluation

Out-of-sample validation was conducted over a 15-day forward horizon, evaluating the structural model against a Naive Baseline Benchmark (1-day lagged return persistence model where predicted return equals yesterday’s actual return).

Performance Metrics Output

• Structural Model RMSE: 0.0241
• Naive Baseline RMSE: 0.0583
• Structural Model R-Squared: 82.08%
• Naive Baseline R-Squared: 4.12%

The Structural Model achieved an out-of-sample R-Squared of 82.08%, significantly outperforming the Naive Baseline. This confirms that equity return variance in MSTR is predominantly explained by spot Bitcoin fluctuations and treasury balance sheet updates rather than simple price momentum.

# ==============================================================================
# SEC XBRL DRIVEN QUANT PIPELINE: MICROSTRATEGY BALANCE SHEET SHOCKS VS BTC
# Author: Selcuk Disci (datageeek.com)
# ==============================================================================

# 1. LOAD REQUIRED LIBRARIES (AUTOMATED PACMAN ENTIRE PIPELINE INGESTION)
# ------------------------------------------------------------------------------
# Check and install pacman package manager if not already available
if (!require("pacman")) install.packages("pacman")
# Ingest core financial, data manipulation, SEC scraping, and modeling libraries
pacman::p_load(secfile, tidyquant, tidyverse, zoo, timetk, tidymodels)

# 2. DEFINE SEC COMPLIANT USER AGENT AND FETCH DATA
# ------------------------------------------------------------------------------
# Set contact email required by SEC EDGAR fair access policy header guidelines
user_agent <- "<your_email_address>"

# Fetch SEC Central Index Key (CIK) identifier for MicroStrategy Inc.
mstr_cik <- get_ciks("MSTR", user_agent = user_agent)

# Retrieve corporate submissions metadata and extract structured XBRL financial facts
mstr_subs <- get_submissions(mstr_cik, user_agent = user_agent)
mstr_facts <- get_data(mstr_subs, user_agent = user_agent)

# 3. EXTRACTION OF BITCOIN HOLDINGS DYNAMICALLY FROM WIDE FORMAT
# ------------------------------------------------------------------------------
# Extract column names from SEC dataset to locate dynamic crypto asset tags
available_columns <- names(mstr_facts)
# Detect target columns matching SEC XBRL taxonomy for digital holdings
target_btc_column <- available_columns[str_detect(available_columns, "DigitalAsset|CryptocurrencyHoldings")]

# Fallback safety net to ensure execution if specific digital asset tags are absent
if(length(target_btc_column) == 0) {
  target_btc_column <- "Assets" 
} else {
  target_btc_column <- target_btc_column[1]
}

# Clean and transform raw Bitcoin holdings time series data
mstr_btc_holdings <- mstr_facts %>%
  select(report_date, !!sym(target_btc_column)) %>% # Dynamically unquote 'target_btc_column' using !!sym() to evaluate the string as a column name
  rename(date = report_date, BTC_Held = !!sym(target_btc_column)) %>%
  mutate(
    date = as.Date(date),
    BTC_Held = as.numeric(BTC_Held) 
  ) %>%
  filter(!is.na(BTC_Held)) %>%
  distinct(date, .keep_all = TRUE) %>%
  arrange(date)

# 4. FETCH AND COMPUTE LOG RETURNS USING TIDYQUANT
# ------------------------------------------------------------------------------
# Define equity and cryptocurrency symbols for market data retrieval
tickers <- c("MSTR", "BTC-USD")

# Download adjusted close daily historical price series and calculate log returns
market_returns <- tq_get(tickers, from = "2020-01-01", get = "stock.prices") %>%
  group_by(symbol) %>%
  tq_mutate(select = adjusted,
            mutate_fun = periodReturn,
            period = "daily",
            type = "log",
            col_rename = "log_return") %>%
  select(date, symbol, log_return, adjusted) %>%
  ungroup()

# Pivot market data into a wide format and standardize asset-specific return/price columns
market_pivoted <- market_returns %>%
  pivot_wider(names_from = symbol, values_from = c(log_return, adjusted)) %>%
  rename(
    MSTR_Log_Return = log_return_MSTR,
    BTC_Log_Return = `log_return_BTC-USD`,
    MSTR_Close = adjusted_MSTR,
    BTC_Close = `adjusted_BTC-USD`
  ) %>%
  mutate(date = as.Date(date))

# 5. ALIGN MIXED FREQUENCY DATA AND COMPUTE SHOCKS
# ------------------------------------------------------------------------------
# Join balance sheet data with market returns and impute missing daily holdings values via forward fill
processed_model_data <- market_pivoted %>%
  left_join(mstr_btc_holdings, by = "date") %>%
  mutate(BTC_Held_Daily = na.locf(BTC_Held, na.rm = FALSE)) %>%
  filter(!is.na(BTC_Held_Daily) & !is.na(MSTR_Log_Return) & !is.na(BTC_Log_Return)) %>%
  mutate(Balance_Sheet_Shock = log(BTC_Held_Daily / lag(BTC_Held_Daily))) %>%
  filter(!is.na(Balance_Sheet_Shock) & is.finite(Balance_Sheet_Shock))

# 6. TIME-BASED DATA SPLITTING VIA TIMETK (STRICT TIME WINDOWS)
# ------------------------------------------------------------------------------
# Partition dataset chronologically to prevent future data leakage during evaluation
data_splits <- time_series_split(
  data       = processed_model_data,
  date_var   = date,
  initial    = "1 year",   
  assess     = "15 days",  
  cumulative = FALSE       
)

# Extract historical training split and testing horizon
train_data <- training(data_splits)
test_data  <- testing(data_splits)

# 7. ESTIMATE LINEAR REGRESSION VIA NATIVE TIDYMODELS COMPONENT
# ------------------------------------------------------------------------------
# Define features profile inside the standard recipe framework
mstr_recipe <- recipe(MSTR_Log_Return ~ date + BTC_Log_Return + Balance_Sheet_Shock, data = train_data) %>%
  update_role(date, new_role = "id")

# Specify parsnip linear regression engine spec
lm_spec <- linear_reg() %>%
  set_engine("lm") %>%
  set_mode("regression")

# Bind graph components into a clean workflow architecture
mstr_workflow <- workflow() %>%
  add_recipe(mstr_recipe) %>%
  add_model(lm_spec)

# Fit model natively on the chronological training split partition
fitted_lm_workflow <- fit(mstr_workflow, data = train_data)

# 8. OUT-OF-SAMPLE PERFORMANCE TESTING WITH NAIVE BENCHMARK USING YARDSTICK
# ------------------------------------------------------------------------------
# Generate clean out-of-sample predictions via standard tidymodels syntax
test_data_predictions <- predict(fitted_lm_workflow, new_data = test_data)

# Construct evaluation frame with actual returns and naive lag benchmark
evaluation_df <- test_data %>%
  select(date, MSTR_Log_Return, MSTR_Close) %>%
  bind_cols(test_data_predictions) %>%
  rename(Predicted_Return = .pred) %>%
  mutate(Naive_Predicted_Return = lag(MSTR_Log_Return, default = first(MSTR_Log_Return)))

# Reshape predictions to comparative long format for unified performance calculation
evaluation_long <- evaluation_df %>%
  select(date, MSTR_Log_Return, Predicted_Return, Naive_Predicted_Return) %>%
  pivot_longer(
    cols = c(Predicted_Return, Naive_Predicted_Return),
    names_to = "model_type",
    values_to = "estimate"
  ) %>%
  rename(truth = MSTR_Log_Return) %>%
  mutate(model_type = if_else(model_type == "Predicted_Return", "Structural_Model", "Naive_Baseline"))

# Compute RMSE and R-Squared accuracy metrics across structural vs naive models
my_financial_metrics   <- metric_set(rmse, rsq)
accuracy_report_tibble <- evaluation_long %>%
  group_by(model_type) %>%
  my_financial_metrics(truth = truth, estimate = estimate) %>%
  ungroup() %>%
  arrange(.metric, model_type)

# Print execution metric table to console
print(accuracy_report_tibble)

# 9. MODERN INTERACTIVE PLOTLY VISUALIZATION (PRICE-BASED DYNAMIC RSI)
# ------------------------------------------------------------------------------
# Load interactive graphics packages required for dynamic reporting
if (!require("pacman")) install.packages("pacman")
pacman::p_load(plotly, scales, glue, ggtext)

# Extract out-of-sample RMSE metric value for confidence band projection
trusted_rmse <- accuracy_report_tibble %>%
  filter(model_type == "Structural_Model" & .metric == "rmse") %>%
  pull(.estimate)

# Extract out-of-sample R-Squared metric value for header annotation
rsq_val <- accuracy_report_tibble %>%
  filter(model_type == "Structural_Model" & .metric == "rsq") %>%
  pull(.estimate)

# Convert predicted log returns back into absolute dollar prices with confidence intervals
df_eval <- evaluation_df %>%
  mutate(MSTR_Yesterday_Close = lag(MSTR_Close, default = first(MSTR_Close))) %>%
  mutate(
    actual  = MSTR_Close,
    pred    = MSTR_Yesterday_Close * exp(Predicted_Return),
    conf_hi = MSTR_Yesterday_Close * exp(Predicted_Return + (2 * trusted_rmse)),
    conf_lo = MSTR_Yesterday_Close * exp(Predicted_Return - (2 * trusted_rmse))
  ) %>%
  filter(date > min(date))

# Build customized hover text data frames for interactive Plotly tooltips
df_plot_actual <- df_eval %>% select(date, actual)  %>% mutate(text_actual = glue("<b>Actual MSTR Price:</b> ${round(actual, 2)}<br><b>Date:</b> {format(date, '%b %d, %Y')}"))
df_plot_pred   <- df_eval %>% select(date, pred)    %>% mutate(text_pred   = glue("<b>Linear AI Pred:</b> ${round(pred, 2)}<br><b>Date:</b> {format(date, '%b %d, %Y')}"))
df_plot_hi     <- df_eval %>% select(date, conf_hi) %>% mutate(text_hi     = glue("<b>Overbought Ceiling:</b> ${round(conf_hi, 2)}<br><b>Date:</b> {format(date, '%b %d, %Y')}"))
df_plot_lo     <- df_eval %>% select(date, conf_lo) %>% mutate(text_lo     = glue("<b>Oversold Floor:</b> ${round(conf_lo, 2)}<br><b>Date:</b> {format(date, '%b %d, %Y')}"))

# Assemble base ggplot layer with confidence bands, actuals, and model projections
p <- ggplot() +
  geom_ribbon(data = df_eval, aes(x = date, ymin = conf_lo, ymax = conf_hi), fill = "#808080", alpha = 0.18) +
  geom_point(data = df_plot_hi, aes(x = date, y = conf_hi, text = text_hi), color = "transparent", alpha = 0, size = 3) +
  geom_point(data = df_plot_lo, aes(x = date, y = conf_lo, text = text_lo), color = "transparent", alpha = 0, size = 3) +
  geom_line(data = df_plot_actual, aes(x = date, y = actual), color = "#2c3e50", linewidth = 1.2) +
  geom_point(data = df_plot_actual, aes(x = date, y = actual, text = text_actual), color = "#2c3e50", size = 2) +
  geom_line(data = df_plot_pred, aes(x = date, y = pred), color = "#e74c3c", linetype = "dashed", linewidth = 1.2) +
  geom_point(data = df_plot_pred, aes(x = date, y = pred, text = text_pred), color = "#e74c3c", size = 2) +
  scale_y_continuous(labels = dollar_format(accuracy = 1)) +
  labs(
    x = "", y = "",
    title = paste0(
      "MicroStrategy (MSTR) <span style = 'color:#2c3e50'>Actual Prices</span> vs ",
      "<span style = 'color:#e74c3c'>Tidymodels Linear AI Forecast</span><br>",
      "<span style='font-size:12px; color:#555555;'>15-Day Trading Horizon | Out-of-Sample R-Squared: ", round(rsq_val * 100, 2), "%</span>"
    )
  ) +
  theme_minimal() +
  theme(plot.title = element_markdown(hjust = 0.5, face = "bold"),
        plot.background = element_rect(fill = "#ffffff", color = NA),
        panel.background = element_rect(fill = "#ffffff", color = NA),
        panel.grid.minor = element_blank())

# Set typography styles for HTML dashboard output
font_family <- list(family = "Roboto Slab, Sans-Serif", size = 16)
label_font  <- list(font = list(family = "Roboto Slab, Sans-Serif", size = 13))

# Convert static ggplot to fully interactive HTML Plotly widget
ggplotly(p, tooltip = "text") %>% 
  style(hoverlabel = label_font) %>% 
  layout(font = font_family) %>% 
  config(displayModeBar = FALSE)

6. Price Transformation & Dynamic Volatility Bands

Log return forecasts are converted back into actionable nominal price predictions using exponential compounding anchored to the previous trading day’s close price. The predicted price equals yesterday’s closing price multiplied by the exponential of the predicted log return.

Empirical Confidence Envelopes (Overbought / Oversold Zones)

Using the empirical Root Mean Squared Error (RMSE) derived from out-of-sample evaluation, dynamic 2-sigma volatility boundaries are constructed around the predicted price path:

• Overbought Ceiling: Calculated as yesterday’s close multiplied by the exponential of the predicted return plus twice the RMSE.
• Oversold Floor: Calculated as yesterday’s close multiplied by the exponential of the predicted return minus twice the RMSE.

Visual Output Analysis

1. Directional Accuracy: The predicted price path (dashed red line) tracks the actual price trajectory (dark solid line) with high fidelity across the 15-day evaluation window.
2. Volatility Envelope Containment: Actual spot equity prices remain fully bounded within the 2-sigma gray confidence envelope, validating the structural model’s error boundary calibration.
3. Price Acceleration Tracking: The sharp upward price repricing between September 14 and September 21 is effectively captured by the structural model due to the immediate integration of spot BTC log returns.

Conclusion

By combining automated SEC EDGAR parsing via secfile package with functional machine learning pipelines in tidymodels, quantitative analysts can build scalable, production-ready asset pricing engines. In the case of MicroStrategy (MSTR), incorporating SEC XBRL balance sheet tracking alongside high-frequency spot crypto returns yields an empirical out-of-sample explanatory power (R-Squared) exceeding 82%.