Skip to content

Instantly share code, notes, and snippets.

View tslumley's full-sized avatar
😐

Thomas Lumley tslumley

😐
View GitHub Profile
@tslumley
tslumley / app.R
Created January 4, 2019 01:13
Shiny app for exploring posterior distributions given surprising data
library(shiny)
# Define UI for application that draws a histogram
ui <- fluidPage(
# Application title
titlePanel("Bayesian Surprise"),
# Sidebar with a slider input for number of bins
@tslumley
tslumley / seven_of_nine.txt
Created November 12, 2023 05:54
Words of nine letters with seven distinct letters
abilities
academics
accessing
accessory
activated
addiction
additions
admission
advantage
aerospace
paired_test<-function(formula, data, subset, paired=NULL, ...){
if (is.null(paired))
return("do the current thing")
## make sure paired is ~id
if (!is.language(paired))
stop("paired must be a formula")
if(!(length(paired)==2 && paired[[1L]]=="~" && length(paired[[2L]])<=2))
@tslumley
tslumley / new_sample_int.R
Last active August 31, 2023 10:28
Modify sample.int with new options for PPS with replacement sampling that have the specified inclusion probabilities
## New sample.int
sample_int<-function(n, size=NULL, replace=FALSE, prob=NULL,
useHash=(n > 1e+07 && !replace && is.null(prob) && (!is.null(size)) && size <= n/2),
method=c("sequential","marginal","poisson")){
if (replace || is.null(prob)){
if (is.null(size)){
size<-n
}
} else{
void tille_incprob(double a[], int *pn, int *plen){
double a_sum=0;
int i,n,len,l,l1;
n=*pn;
len=*plen;
for(i=0;i<len; i++){
a_sum+=a[i];
}
@tslumley
tslumley / pairwise.R
Last active March 31, 2023 00:28
pairwise lmm with kinship
Documentation
X: matrix of predictors (including the intercept)
y: matrix of outcomes
kinship: the kinship matrix (with 1s on the diagonal, not 2s)
pwt_mat: matrix of pairwise sampling weights (ie, reciprocal of pairwise sampling probabilities)
Output:
1. genetic variance divided by residual variance
2. residual variance
@tslumley
tslumley / brant_test.R
Created February 11, 2023 04:31
Based on Brant's test for the proportional odds assumption
olr_brant_test<-function(formula, design,test=c("brant-original","omnidirectional-Wald")){
test<-match.arg(test)
m1<-svyolr(formula, design = design)
K<-length(m1$lev)
P<-length(m1$coef)
get_infl<-function(k,formula,design){
y<-formula[[2]]
formula[[2]]<-bquote(I(as.numeric(.(y))>.(k)))
mk<-svyglm(formula, design, family=quasibinomial, influence=TRUE)
@tslumley
tslumley / busdata.txt
Last active November 2, 2022 22:17
Buses
mo day cancels
11 2 1806
11 1 1846
10 31 1887
10 30 1276
10 29 914
10 28 1828
10 27 1672
10 26 1686
10 25 1845
@tslumley
tslumley / pdl.R
Last active October 13, 2022 06:45
Polynomial distributed lag glms
lagged<-function(x,lag=1){
if (lag==0) return(x)
n<-length(x)
c(rep(NA,lag),x[-( (n-lag+1):n)])
}
pdlweights<-function(lag,degree,tiedown=c(F,F)){
if (any(tiedown)) stop("Tiedown not working")
contr.poly(lag,contrasts=F)[,tiedown[1]+(1:degree)]
@tslumley
tslumley / README
Created July 23, 2022 08:37
vctrs class for multiple-response objects
Based on the `rimu` package, the start of an implementation of 'mr' objects using Hadley's `vctrs` package, so they can be used in tibbles.