Ben Bolker

@bbolker.bsky.social

2401 Followers 631 Following

Ecology, evolution, epidemiology, statistics (mixed models). McMaster University, Hamilton, Ontario https://math.mcmaster.ca/bolker

set.seed(101) library(glmmTMB) library(lme4) dd <- expand.grid(f = factor(rep(LETTERS[1:4])), g = factor(1:25), rep = 1:10) betavec <- rnorm(4) sdvec <- rnorm(4) corvec <- rnorm(choose(4,2)) dd$y <- simulate_new(~f + us(f|g), family = gaussian, newparams = list(beta = betavec, betadisp = -2, theta = c(sdvec, corvec)), newdata = dd)[[1]] m1 <- glmmTMB(y ~ f + us(f|g), data = dd, family = gaussian) ## FIXME: warn if additional args are specified? m2 <- glmmTMB(y ~ f + homcs(0+f|g), data = dd, family = gaussian) m3 <- lmer(y ~ f + cs(0+f|g, hom = TRUE), data = dd, REML = FALSE) stopifnot(all.equal(logLik(m2), logLik(m3))) m4 <- glmmTMB(y ~ f + (1|g/f), data = dd, family = gaussian) m5 <- lmer(y ~ f + (1|g/f), data = dd, REML = FALSE) stopifnot(all.equal(logLik(m2), logLik(m4))) stopifnot(all.equal(logLik(m5), logLik(m4)))
If worn correctly, the small box respirator, introduced by the British late in 1916 (though not fully issued until January 1917),41 provided reliable protection against chlorine and phosgene.42 Largely for this reason, in July 1917 the Germans introduced mustard gas, a vesicant that burned the skin. Known as ‘yellow cross’ from the shell markings, the toxin immediately produced a flood of admissions: 160,970 soldiers presented at CCSs between 21 July and 23 November 1917, of whom 1,859 died.43 It was estimated that 75 per cent had been exposed to mustard gas.44 Deaths initially occurred because the toxin had delayed effects.45 The smell of the gas was ‘not very noticeable and… the immediate effect on troops exposed to it’, Brigadier K. Wigram reported, was ‘only a slight irritation of the nose and throat’, so that most troops had been slow to realize that they had been poisoned.46 Mustard gas achieved its greatest effect in the months immediately after its introduction because ‘the novelty of the condition, the multiplicity of the symptoms, and the entire absence of knowledge as to possible after-effects naturally led to the condition being overestimated’.47 Fear inspired by mustard gas was heightened by its capacity to cause loss of sight. Although most blinded servicemen recovered, the acute photophobia, conjunctivitis, and oedema of the eyelids forced soldiers to close their eyes, ‘so much so in fact that when some of the milder cases were evacuated each man had to be led like a blind man by an orderly to the ambulance car’.48 Such images struck a cultural chord and remain among the iconic symbols of the First World War.49
To put things less delicately, Chernoff faces are an exam- ple of an idea that is widely acknowledged to be bad, that is used by essentially nobody in the real world, and yet is nonetheless widely circulated in academia. Contrast this with rainbow colormaps and 3D pie charts, which have also been generally complained about, but are still seen in the wild with some frequency.
side-by-side plots: (A) a traditional linegraph of bitcoin prices (a jagged, random-walk looking line); (B) automatically restyled as a painting, with white mountains matching the trajectory of the line, with a blue sky above (with stars added) and a gray-green foreground. Caption: "Figure 7: Using a neural net to automatically restyle a time series, from Moritz Klack [14]. This would be a useful way of generating custom Ross-Chernoff plots, and also allows me to include “deep learning” in the author keywords"
Kobo's "young adult books on sale" showing thumbnails for _The Chronicles of Narnia_, _Brave New World_, _Mission Improbable_, and _The Stranger_
Rather than automatically removing predictors, focus on whether multicollinearity prevents you from answering your specific research question, and whether the instability in coefficient estimates is acceptable given your goals.
A line plot of four (inverse) link functions [inverse-](cauchit, cloglog, logit, probit) for x from -5 to 5. Probability is represented on the logit scale, so the (inverse)logit is linear by definition. Cauchit is flat and symmetric. Probit has a slightly different scale and is also curvy in the tails. Cloglog increases very rapidly above x=0, and is capped at 1.0 for values > 3. Code: library(tidyverse); theme_set(theme_bw()) xvec <- seq(-5, 5, length = 101) links <- c("logit", "probit", "cauchit", "cloglog") dd <- setNames(links, links) |> purrr::map_dfr(~ data.frame(x = xvec, y = make.link(.)$linkinv(xvec)), .id = "link") brkvec <- c(1e-5, 0.5, 0.9, 0.999999) ggplot(dd, aes(x, y, color = link)) + geom_line() + scale_colour_brewer(palette = "Dark2") + scale_y_continuous(trans = "logit", breaks = brkvec, name = "probability (logit scale)") ggsave("links.png")
INCUBATION PERIOD: 2 to 4 weeks (range from few days to 2 months) for HFRS; and 14-17 days for HPS Footnote 3Footnote 4. COMMUNICABILITY: Person-to-person transmission is very rare, but has been observed for Andes virus during an outbreak in Southern Argentina Footnote 1Footnote 7.
A snippet of the Muncie Evening Press (p. 4) from Weds Feb 17 1932. We can see the words "... mostly on experience and ... experience usually comes from ... judgment"
Panel 1 Caption: “SO THE THREE HAPPY FRIENDS WERE JUST STANDING THERE, MINDING THEIR OWN BUSINESS…” Three simple figures stand side by side. Their torsos are labeled, from left to right: “eˣ,” “eʸ,” and “yˣ.” Panel 2 Caption: “…WHEN SUDDENLY… OH NO… HERE COMES THE ∂/∂x MONSTER!” A large dinosaur‑like creature labeled “∂/∂x” approaches from the left toward the three figures, who look alarmed. Panel 3 Caption: “OMG! HE’S EATING MISTER eˣ!” The derivative monster bites and swallows the figure labeled “eˣ.” The figures labeled “eʸ” and “yˣ” watch with horrified expressions. Panel 4 Caption: “BUT WAIT… IT’S OK. HE’S JUST GOING TO POOP HIM RIGHT BACK OUT UNHARMED, THANK GOD.” The monster crouches and poops out the figure labeled “eˣ”, who bounces once and lands on the floor, smiling. Panel 5 Caption: “BUT, OH NO, HERE HE COMES AGAIN! RUN! RUN AWAY!… DAMN, IT’S TOO LATE. HE’S GOT MR eˣ!” The monster lunges back toward the “eˣ” figure and grabs it. Mr. "yˣ" looks on in horror. Panel 6 Caption: “I CAN’T BEAR TO WATCH THIS. HE’S MAULING MR eˣ TO DEATH! SWEET FANCY MOSES!” The monster lifts and chews the “eˣ” figure. Mr. "yˣ" watches. Panel 7 Caption: “THERE’S eˣ BLOOD AND GUTS SPEWING ALL OVER THE PLACE! OH, THE HUMANITY!!!” Exaggerated gore—detached figure parts and splashes—erupts as the monster continues. Panel 8 Caption: “RUN, MR yˣ, RUN! THERE’S NO TELLIN’ WHAT THAT ∂/∂x BASTARD IS GONNA DO TO YOUR SORRY ASS…” The “yˣ” figure runs to the right. A speech bubble from off‑panel says: “OK, BILLY, THAT’S ENOUGH. SIT DOWN.” Panel 9 The scene is revealed to be a classroom. Billy is working on a problem on the board while his teacher (female presenting) looks on. He says: "BUT I'M NOT FINISHED YET." She says: "YES YOU ARE". Chalkboard text (visible): “∂/∂x ( eˣ + eʸ + yˣ )” “= eˣ + 0 + …”
A State Department official told Fox News Digital that the new designs will be available for "any American citizen" who applies for a passport when the rollout happens and will continue for as long as there is availability. The passports will only be available at the Washington Passport Agency, the official detailed.
same picture as before, but with overseas territories circled in red (two in the Pacific, at left edge; four in South Atlantic; two in north Atlantic [Bermuda, Gibraltar])
a map of the UK that includes a bunch of overseas possession in the Atlantic and Pacific, so Great Britain is a very small blob in the upper right hand corner.
A page of complicated sheet music that is more or less gibberish
image showing dates of past NSERC Discovery Grant announcements, with three GAMs (cubic spline, Gaussian process, thin plate). x-axis years 2013-2026 (data up to 2025), y-axis announcement date (March 15 - April 22). 2013 date is around March 22, curve peaks in 202 around April 17. Prediction for GP in 2026 is March 30 (2 SE = 22 March - 7 April). uncertainty estimates for cubic spline, thin plate spline are wider (March 15 - April 15)
a picture of a white cat peering from behind a computer screen. The screen is showing part of a Firefox window; only the left half of the cat's face is visible.
A Shiny app showing a subset of lines selected from a prior sample
A screenshot from the Merlin app showing Bald Eagle, Eastern Meadowlark, Killdeer, and Dark-eyed Junco
## fit multi-type model (don't scale count var) form_MMLM2 <- cbind(femur_s, tibia_s, tarsus_s, SCT) ~ trait:(genotype*temp) - 1 m2 <- MCMCglmm(form_MMLM2, random=~ us(trait):line, rcov=~ us(trait):units, data= as.data.frame(dll_data), family = rep(c("gaussian", "poisson"), times = c(3, 1)), nitt= 10000, burnin= 2000, thin=10)
badge from web site: "CRAN synchronization failing"
The Sakura Cherry Blossoms in High Park are firmly closed in brown buds. It is still too early to make a definitive prediction on when peak bloom will occur. Note: Sakura Cherry Blossoms in High Park typically reach peak bloom anytime between late April to early May.
For works first published prior to 1978, the term will vary depending on several factors. To determine the length of copyright protection for a particular work, consult chapter 3 of the Copyright Act (title 17 of the United States Code). More information on the term of copyright can be found in Circular 15a, Duration of Copyright, and Circular 1, Copyright Basics.
The Bayesian results imply much higher risk of early collapse than maximum likelihood methods. This difference is due to large probabilities of early collapse for certain parameter values that are plausible in light of the data. Because of simplifying assumptions, these results are not directly applicable to assessment. Nevertheless they imply that maximum likelihood and similar methods based upon point parameter estimates will grossly underestimate the risk of early collapse.
"Embryologist E. E. Just complained that genetics and selection could explain why populations of flies had more or fewer bristles on their backs, but it couldn't explain how a fly constructed its back in the first place (Harrison 1937: 372; Gilbert et al 1996: 361)
results of `apropos("R.?[Vv]ersion", ignore.case=FALSE)`: c("getRversion", "R.version", "R.Version", "R.version.string")
Circle 4 Over-Vectorizing We skirted past Plutus, the fierce wolf with a swollen face, down into the fourth Circle. Here we found the lustful. It is a good thing to want to vectorize when there is no effective way to do so. It is a bad thing to attempt it anyway. A common reflex is to use a function in the apply family. This is not vectorization, it is loop-hiding. The apply function has a for loop in its definition. The lapply function buries the loop, but execution times tend to be roughly equal to an explicit for loop. (Confusion over this is understandable, as there is a significant difference in execution speed with at least some versions of S+.) Table 4.1 summarizes the uses of the apply family of functions. Base your decision of using an apply function on Uwe’s Maxim (page 20). The issue is of human time rather than silicon chip time. Human time can be wasted by taking longer to write the code, and (often much more importantly) by taking more time to understand subsequently what it does
Image from a 14th c. Italian manuscript of the Inferno in the Bodleian library at Oxford https://digital.bodleian.ox.ac.uk/objects/ab35e336-a471-4cf0-a9a7-592dbb8695d8/surfaces/cba2ff2c-ecd7-44e3-907b-1bd6d0376daf/ : Dante and Virgil are at the right of the image, facing left towards the minotaur who is running down a hill. In this version the head and torso are a man's, the body is the bull's.
Left panel, a salt shaker containing multiple pepper packets, labeled 'x'. Middle, shaker containing one packet, labeled 'x[0:1]'. Right, top, one pepper packet, labeled 'x[0]'. Right, bottom, pile of pepper grains, labeled 'x[0][0]'
What’s happening? dd[["x"]] and dd$x return the vector integer64, so sum() dispatches to bit64’s method for the Summary group generics (sum.integer64), which knows how to add 64‑bit integers correctly. Result: 55 ✅ dd["x"] returns a one-column data.frame (i.e., a list). When you call sum(dd["x"]), you’re now calling sum() on a data.frame, not on an integer64 vector. For data.frames, base R uses the Summary group generic method Summary.data.frame (in base). Internally, this method coerces the data frame to something vector-like (e.g., via as.matrix() or unlist() pathways) before applying the sum generic. In doing so, it strips the integer64 class, and the remaining underlying storage is misinterpreted. integer64 (from bit64) is not a native R type; it’s a classed vector with an underlying storage that does not map 1:1 to base numeric types. When the class is stripped by the data.frame method, the bytes get reinterpreted (e.g., as doubles) instead of being properly converted. That’s why you see a tiny, nonsensical denormal-like value such as 2.717361e-322 — it’s the result of interpreting the 64‑bit integer bits as an IEEE-754 double, not the sum of actual numbers.
> library(bit64) > dd <- data.frame(x = as.integer64(1:10), y = as.integer(1:10)) > sum(dd$x) integer64 [1] 55 > sum(dd["x"]) [1] 2.717361e-322 > sum(dd[["x"]]) integer64 [1] 55 > sum(dd$y) [1] 55 > sum(dd["y"]) [1] 55 > sum(dd[["y"]]) [1] 55 >
grepv("^c", rownames(installed.packages())) |> Filter(f = \(x) x %in% qdap::GradyAugmented)
R code library(microbenchmark); library(ggplot2) set.seed(1001) n <- runif(1e6, 0, 2) m1 <- microbenchmark( pmax(n,1), pmax.int(n,1), {n[n<1] <- 1} ) autoplot(m1)
A line graph with sample size per group on the x axis, proportion of simulations on the y-axis (empirical estimate of probability). Lines disaggregate outcomes: "clearly large and positive" starts at about 0.08 for n = 5, grows to near 1.0 for n>100. "Clearly positive, maybe large" peaks around n=10. Not both (large and negative) declines from 0.3 to nearly 0 at n=20. Other categories are very low prob in this case.
A flow chart from Schruben 1980 showing (top) 'exogenous inputs'. These flow to 'real system' (left) and 'simulated system' (right). Both the real and simulated systems flow to a report; the two reports then flow into a 'shuffle' step, which then goes to 'manager', and finally (bottom) splits to two boxes labeled 'real' and 'bogus'.
screenshot from the Ontario respiratory virus dashboard showing % positivity [COVID-19, influenza, RSV]; cases [COVID-19, influenza]; outbreaks [all three]; hospital admissions [all 3]. Focusing on influenza, the top three indicators (% positivity, cases, outbreaks) all appear to have peaked recently
a 3D surface plot showing (x) prey size, (y) prey density, (z) proportion predated
Cytokine-mediated CAR T therapy resistance in AML Anand S. Bhagwat, Leonel Torres, Olga Shestova, Maksim Shestov, Patrick W. Mellors, Hansen R. Fisher, Saamia N. Farooki, Benjamin F. Frost, Michael R. Loken, Avery L. Gaymon, Diane Frazee, Walter Rogal, Noelle Frey, Elizabeth O. Hexner, Selina M. Luger, Alison W. Loren, Mary Ellen Martin, Shannon R. McCurdy, Alexander E. Perl, Edward A. Stadtmauer, Jennifer L. Brogdon, Joseph A. Fraietta, Wei-Ting Hwang, Don L. Siegel, Gabriela Plesa, Richard Aplenc, David L. Porter, Carl H. June & Saar I. Gill
library(ggplot2) data("sleepstudy", package = "lme4") ggplot(sleepstudy, aes(Days, Reaction)) + geom_point() + stat_summary(geom="point", fun = mean, color = "red") + stat_summary(geom="line", fun = mean, color = "red") ggsave("linespoints.png")
Curve showing relative sunset times for different locations in southern Canada/northeastern US. Parabolic curves with minima at December 7-8, y-axis shows an increase of 0.4 minutes in the sunset time as go back to December 4 or forward to December 12
a list of package maintained by John Fox, from https://cran.r-project.org/web/checks/check_summary_by_maintainer.html
I’m working on a manuscript that mentions a fair number of standard distributions and their properties, and I’ve been wondering if there’s an equally standard reference. For context for my question, while Wikipedia often has what I need if I want to appeal to some function and its properties, the Digital Library of Mathematical Functions (https://dlmf.nist.gov/) provides a peer-reviewed standard alternative. Is there someplace similar where I could point for properties of a Gumbel, Dirichlet-multinomial distribution, etc.
Statisticians are increasingly posed with thought-provoking and even paradoxical questions, challenging our qualifications for entering the statistical paradises created by Big Data. By developing measures for data quality, this article suggests a framework to address such a question: “Which one should I trust more: a 1% survey with 60% response rate or a self-reported administrative dataset covering 80% of the population?” A 5-element Eulerformula-like identity shows that for any dataset of size n, probabilistic or not, the difference between the sample average Xn and the population average XN is the product of three terms: (1) a data quality measure, ρR,X, the correlation between Xj and the response/recording indicator Rj ; (2) a data quantity measure, √(N − n)/n, where N is the population size; and (3) a problem difficulty measure, σX, the standard deviation of X. This decomposition provides multiple insights: (I) Probabilistic sampling ensures high data quality by controlling ρR,X at the level of N −1/2; (II) When we lose this control, the impact of N is no longer canceled by ρR,X, leading to a Law of Large Populations (LLP), that is, our estimation error, relative to the benchmarking rate 1/√n, increases with √N; and (III) the “bigness” of such Big Data (for population inferences) should be measured by the relative size f = n/N , not the absolute size n; (IV) When combining data sources for population inferences, those relatively tiny but higher quality ones should be given far more weights than suggested by their sizes.
Note on data-based priors Because the scaling is based on the scales of the predictors (and possibly the outcome) these are technically data-dependent priors. However, since these priors are quite wide (and in most cases rather conservative), the amount of information used is weak and mainly takes into account the order of magnitude of the variables. This enables rstanarm to offer defaults that are reasonable for many models.
Time series of CAD from 1950 to 2025. Value fluctuated between $0.90 and $1.05 USD from 1950-1975, $0.60 to $0.90 1975-2005, briefly exceeded par around 2010, now back down to about $0.75 USD
Keron Thomas (c. 1975 – July 22, 2013)[2] was primarily known for posing as New York City Subway motorman "Regoberto Sabio" on May 8, 1993, and operating an A train in revenue service for over three hours, when he was sixteen years old.
plot of scraped data from Edwards et al 2021. On a logarithmic y scale, regression slopes are very similar. Red points are superspreaders (nominally > 156), black are 'normal' (nominally < 156)
 ‘What would you do,’ I had asked, ‘if drawing a Latin square at random for an experiment, you happened to draw a Knut-Vik (Knight’s move) square?’ Sir Ronald said he thought he would draw again and that, ideally, a theory explicitly excluding regular squares should be developed.” (L. J. Savage to R. A. Fisher, 1952, p. 88, ref. 21).
the title block of an obviously AI-generated paper on Authorea (although maybe only AI-organized?): Certainly! Apologies for the previous omissions. Below is the complete LaTeX document that includes all the requested sections, arguments, code snippets, and proofs, organized logically into a single cohesive document. “‘latex Diagnosis and Treatment of Foot Fusarium Infection Confirmed by Next-Generation Sequencing Technology: A Case Report and Literature Review dermatology +7 Yulin Wang, Danyue He, Qi Zhao, Jie Hu, Ruyi Wang, Yiyang Hu, Xiangning Li, Jia Cai, Lili Tang, Xing Fan
Figure 2: Four graphs showing patterns of causal effects, each with average effect of 0.1, but varying in different ways as a function of a pre-treatment predictor: (e) linear interaction, (f ) no effect then steady increase, (g) plateau, (h) intermediate zone with large effects
a map of the incidence of favism (G6PD deficiency), showing concentrations in sub-Saharan Africa and the middle east, some in India, SE asia, and Brazil
library(vroom) system.time({ dd <- rnorm(15e6*10) |> matrix(nrow = 15e6) |> as.data.frame() vroom_write(dd, "tmp.csv") dd2 <- vroom("tmp.csv") } ) ## user system elapsed ## 34.995 5.127 8.327
Lucy (L) and Charlie Brown (R) are standing behind a wall. (Panel 1) CB: "I hate being a nothing! I refuse to go through the rest of my life as a zero!" (2) L: "What would you like to be, Charlie Brown, a five? Or how about a twenty-six? Or a par seventy-two?" (3) L: "I know what you could be, Charlie Brown ... a square root!" (4) L: "I think you'd make a great square root, Charlie Brown ..." CB: "I can't stand it!"
A screenshot of samples of digits of the Euler-Mascheroni constant (go to the URL for text-based version)
An image of beak width vs beak depth for a variety of penguin species.
Based on my search, I cannot find credible medical experts who currently contend that the MMR vaccine should not be given routinely to children. Here's what I found: Lack of Current Expert Opposition: The medical literature and recent sources don't identify specific credentialed medical experts who oppose routine MMR vaccination based on scientific evidence. Studies show that 93% of surveyed physicians agree with current official vaccination recommendations and would apply them to their own children Measles 2025 | New England Journal of Medicine Historical Context: The most prominent medical figure associated with MMR vaccine opposition was Andrew Wakefield, but as mentioned earlier, his 1998 study was retracted for fraud, and he lost his medical license. His work has been thoroughly debunked. Current Political vs. Medical Opposition: While there are concerns about potential policy changes under Health Secretary Robert F. Kennedy Jr., who has historically been vaccine-skeptical medRxiv, this represents political rather than credentialed medical expert opinion. Notably, RFK Jr. recently endorsed the MMR vaccine, which "stoked fury among his supporters" The Evidence on Vaccines and Autism | Johns Hopkins | Bloomberg School of Public Health in the anti-vaccine community. The Distinction: While there are individuals who call themselves doctors or health practitioners who oppose vaccines, finding credentialed medical experts (licensed physicians with relevant training) who oppose routine MMR vaccination based on scientific evidence is extremely difficult because the scientific consensus is so strong in favor of the vaccine. The overwhelming medical and scientific consensus continues to support routine MMR vaccination for children, with major medical organizations, pediatricians, and public health experts maintaining this position based on extensive safety and efficacy data.
Another things that jumps out is that maybe I should drop the early data. First off, there’s lots of reasons to ignore them anyways, since testing infrastructure was so different then, and second off, they make the range of the graph go really low, below minus nine “nels.” So, what is a nel? You’ve probably guessed, it’s called that because it is one Natural Log unit, (one NL; one “nel”). A nel is a big number: 1 nel is over 2.5 times my reference rate, and 2 nels is over 7 times my reference rate. Once I drop the very early data, I might see that a “nel”, like a “bel”, might be an awkwardly sized unit for the kinds of changes I’m most interested in, so we’ll look at “centinels” instead. This is the same idea as changing from a proportion to a percent.
a flow chart showing things you can use LLMs for (learning, getting good marks, getting a job, getting into graduate school). Some of these flow downstream (learning stuff gets good marks; good marks get good jobs/grad school entrance; grad school gets good jobs), all of them (except LLMs) may also lead to happiness. There is a red X on the direct path from LLMs to good marks; this is what we want to forbid.
VSCode installed-extensions tab showing Container Tools, Dev Containers, R, R Syntax
flow chart for assessment of exposure to bats. post-exposure management is only recommended if direct exposure is known to have occurred, or if the bat was in the same room as a person *and* there was evidence of direct contact (i.e. scratch, crying/upset child) with a person who may not be able to give a reliable history (child or cognitively impaired adult)