bbolker

#https://dynamicecology.wordpress.com/2021/03/15/how-long-do-institutional-investigations-into-accusations-of-serious-scientific-misconduct-typically-take-heres-some-data/
## https://secretariat.mcmaster.ca/app/uploads/Research-Integrity-Policy.pdf
## https://www-nature-com.libaccess.lib.mcmaster.ca/articles/d41586-020-00287-y

library(tidyverse)
dd <- read.table(header=TRUE, text="
subject duration year
Fuji 24 2010
Fuji 3 2012
Boldt 21 2010

bbolker

library(lme4)
set.seed(101)
nn <- 100
ng <- 10
dd <- data.frame(x = rnorm(nn),
                 f = factor(rep(1:(nn/ng), ng)),
                 g = factor(rep(1:ng, each = ng)))
dd$y1 <- simulate( ~ x + (1|f),
                 newdata = dd,
                 newparams = list(beta = c(1,1),

bbolker

library(tidyverse)
library(png)
## https://stackoverflow.com/questions/9917049/inserting-an-image-to-ggplot2
## dates are approximately when university research offices notify recipients
##  (before ResearchNet posting, dates may vary by university)
this_year <- 2024
(df <- read.table(text="
2025 April 1
2024 March 27
2023 March 29

bbolker

N <- 100
set.seed(101)
dd <- data.frame(
    y = rpois(N, exp(rnorm(N))),
    f = factor(seq(N))
)
dd <- dd[order(dd$y), ]
library(lme4)
m <- glmer(y ~ 1 + (1|f), family = poisson, data = dd)

bbolker

## https://forward.com/fast-forward/389673/finally-mathematical-proof-that-dreidel-is-a-terrible-game/
## Feinerman, Robert. “An Ancient Unfair Game.” American Mathematical Monthly 83, no. 8 (1976). https://www.jstor.org/stable/2319887.

library(emdbook)
library(ggplot2)
library(rayshader)
library(colorspace)

pvec <- 2:5
nvec <- 2:20

bbolker

library(lme4)
library(sjPlot)
set.seed(101)
n <- 1e3 ## obs per group
ns <- 4  ## groups
dd <- data.frame(s = rep(1:ns, each= n),
                 f = factor(sample(1:20, size = n*ns, replace = TRUE)),
                 x = rnorm(n*ns))
sims <- lapply(1:ns,
               function(i) simulate(~x + (1|f),

bbolker

## hinges as computed by stat_boxplot
my_hinges <- function(y, type = 7, coef = 1.5) {
    qs <- as.numeric(stats::quantile(y, c(0.25, 0.75), type = type))
    iqr <- diff(qs)
    outliers <- y < qs[1]-coef*iqr | y > qs[2] + coef*iqr
    hinges <- range(c(qs, y[!outliers]))
    return(hinges)
}

f_hinges <- function(x) {

bbolker

library(lme4)
library(glmmTMB)
library(glmulti)
set.seed(101)

## updated 3 March 2024; take random effects as formulas so that e.g. `ar1()` terms work

## FIXME: better data so we don't get convergence weirdness etc.?
## A random vector of count data
vy1 <- round(runif(100, min=1,max=20)*round(runif(100,min=1,max=20)))

bbolker

## Comparison of bias-correction methods

## n.b. both "logit-normal" and "logistic-normal" are used
##  for logit(Y) ~ Gaussian(mu, sigma)
library(logitnorm)
library(emmeans)

##
eta <- -2
ranef_sd <- 1

bbolker

library(brms)
library(lme4)
system.time(b1 <- brm(Reaction ~ Days + (Days|Subject), sleepstudy)) ## 39 sec

scode <- make_stancode(Reaction ~ Days + (Days|Subject), data = sleepstudy)
sdata <- make_standata(Reaction ~ Days + (Days|Subject), data = sleepstudy)

library(microbenchmark)
m1 <- microbenchmark(
    lmer = lmer(Reaction ~ Days + (Days|Subject), sleepstudy),