"compare_lineups" <- function(l1.o, l1.d, l2.o, l2.d, output=TRUE) { # calculate mean net pts per poss l1.mean <- mean(l1.o)-mean(l1.d) l2.mean <- mean(l2.o)-mean(l2.d) # calculate std error of mean net pts per poss # get sample variances (in other words, se^2) l1.o.var <- var(l1.o)/NROW(l1.o-1) l1.d.var <- var(l1.d)/NROW(l1.d-1) l2.o.var <- var(l2.o)/NROW(l2.o-1) l2.d.var <- var(l2.d)/NROW(l2.d-1) # get se for mean net pts per poss l1.se <- sqrt( l1.o.var + l1.d.var ) l2.se <- sqrt( l2.o.var + l2.d.var ) net.se <- sqrt( l1.se^2 + l2.se^2 ) # calculate z statistics l1.z <- (l1.mean-l2.mean)/net.se l2.z <- (l2.mean-l1.mean)/net.se # estimate probability lineup1 is better than lineup2, and vice-versa l1.prb <- pnorm(l1.z) l2.prb <- pnorm(l2.z) #l1.prb <- 1-pnorm(0, mean=l1.mean-l2.mean, sd=net.se) #l2.prb <- 1-pnorm(0, mean=l2.mean-l1.mean, sd=net.se) # calculate 95% confidence intervals l1.ci <- qnorm(c(0.025,0.975), mean=l1.mean-l2.mean, sd=net.se) l2.ci <- qnorm(c(0.025,0.975), mean=l2.mean-l1.mean, sd=net.se) if (output) { # show results for L1-L2 cat("For Lineup 1 - Lineup 2:\n",sep="") cat("--> Mean: ",sprintf("%.1f",(l1.mean-l2.mean)*100),"\n",sep="") cat("--> Std Err: ",sprintf("%.1f",net.se*100),"\n",sep="") cat("--> 95% CI: (",sprintf("%.1f",l1.ci[1]*100),", ",sprintf("%.1f",l1.ci[2]*100),")\n",sep="") cat("--> Z-score: ",sprintf("%.2f",l1.z),"\n",sep="") cat("--> Pr(L1 > L2): ",sprintf("%.4f",l1.prb),"\n",sep="") } return( list(diff=(l1.mean-l2.mean)*100, se=net.se*100, l1.z=l1.z, l2.z=l2.z, l1.prb=l1.prb, l2.prb=l2.prb, l1.ci=l1.ci*100, l2.ci=l2.ci*100) ) } # Example from 06-07 #res <- compare_lineups(c(rep(0,39+45),rep(1,2+5),rep(2,29+37),rep(3,2+5)), c(rep(0,52+42),rep(1,3+2),rep(2,31+26),rep(3,6+1)), c(rep(0,87+30),rep(1,15+4),rep(2,66+23),rep(3,9+2)), c(rep(0,29+80),rep(1,3+4),rep(2,18+70),rep(3,9+18)))