To cite the hyper2 package in publications, please use Hankin (2017).
Here I analyse a dataset on NZ university rankings, taken from the
Times Higher Educational Supplement. First load the data:
uni_table <- read.table("universities.txt",header=TRUE)
universities_table <- uni_table
head(universities_table,12)
## uni year Teaching Research Citations Industry International
## 1 UoA 2012 28.3 35.0 55.7 76.6 92.9
## 2 Otago 2012 24.5 28.0 48.0 27.8 88.8
## 3 AUT 2012 DNP DNP DNP DNP DNP
## 4 Canterbury 2012 17.3 24.3 30.5 26.6 76.3
## 5 Lincoln 2012 DNP DNP DNP DNP DNP
## 6 Massey 2012 17.4 17.2 24.5 DNP 77.8
## 7 Victoria 2012 16.5 21.5 51.1 25.8 89.0
## 8 Waikato 2012 13.3 13.9 47.8 24.5 87.0
## 9 UoA 2013 34.0 40.9 64.2 76.6 88.0
## 10 Otago 2013 28.8 30.9 54.8 37.2 85.4
## 11 AUT 2013 DNP DNP DNP DNP DNP
## 12 Canterbury 2013 19.9 24.0 48.4 49.8 88.0
Each row corresponds to a ranking in the THE ranking system. Thus the first row refers to The University of Auckland 2012, and we see that they scored 28.3 for their teaching, 25.0 for their research, and so on; rows 3 and 4 show that AUT and Lincoln were not given rankings that year. Rows 1-8 give results for 2012, rows 9-16 give 2013, and so on to 2020.
We can coerce each year to a hyper2 likelihood function by
converting to an order table with order(), thence to a support function:
f <- function(a,y){
a <- subset(a,a$year==y)
o <- function(v){ # 'v' is a character vector
good <- v != "DNP"
v[good] <- rank(-as.numeric(v[good]),ties.method="last")
return(v)
}
out <- data.frame(o(a[,3]),o(a[,4]),o(a[,5]),o(a[,6]),o(a[,7]))
colnames(out) <- colnames(a)[3:7]
rownames(out) <- a[,1]
return(out)
}
f(universities_table,2012)
## Teaching Research Citations Industry International
## UoA 1 1 1 1 1
## Otago 2 2 3 2 3
## AUT DNP DNP DNP DNP DNP
## Canterbury 4 3 5 3 6
## Lincoln DNP DNP DNP DNP DNP
## Massey 3 5 6 DNP 5
## Victoria 5 4 2 4 2
## Waikato 6 6 4 5 4
Now:
universities <- hyper2(pnames=universities_table[1:8,1])
for(i in 2013:2020){universities %<>% add(ordertable2supp(f(universities_table,i)))}
universities
## log(AUT^13 * (AUT + Canterbury)^-1 * (AUT + Canterbury + Lincoln + Massey)^-5 * (AUT +
## Canterbury + Lincoln + Massey + Otago)^-2 * (AUT + Canterbury + Lincoln + Massey + Otago
## + UoA + Victoria + Waikato)^-40 * (AUT + Canterbury + Lincoln + Massey + Otago +
## Victoria)^-3 * (AUT + Canterbury + Lincoln + Massey + Otago + Victoria + Waikato)^-27 *
## (AUT + Canterbury + Lincoln + Massey + Otago + Waikato)^-7 * (AUT + Canterbury + Lincoln
## + Massey + UoA + Victoria + Waikato)^-3 * (AUT + Canterbury + Lincoln + Massey +
## Victoria)^-3 * (AUT + Canterbury + Lincoln + Massey + Victoria + Waikato)^-17 * (AUT +
## Canterbury + Lincoln + Massey + Waikato)^-15 * (AUT + Canterbury + Lincoln + Otago +
## Victoria + Waikato)^-1 * (AUT + Canterbury + Massey + Otago)^-2 * (AUT + Canterbury +
## Massey + Otago + Waikato)^-3 * (AUT + Canterbury + Massey + Victoria + Waikato)^-3 * (AUT
## + Canterbury + Massey + Waikato)^-2 * (AUT + Canterbury + Otago)^-2 * (AUT + Lincoln)^-6
## * (AUT + Lincoln + Massey)^-10 * (AUT + Lincoln + Massey + Otago)^-2 * (AUT + Lincoln +
## Massey + Otago + UoA + Victoria + Waikato)^-3 * (AUT + Lincoln + Massey + Otago +
## Victoria)^-1 * (AUT + Lincoln + Massey + Otago + Victoria + Waikato)^-5 * (AUT + Lincoln
## + Massey + Otago + Waikato)^-3 * (AUT + Lincoln + Massey + Victoria)^-2 * (AUT + Lincoln
## + Massey + Victoria + Waikato)^-2 * (AUT + Lincoln + Massey + Waikato)^-14 * (AUT +
## Lincoln + Otago + Victoria + Waikato)^-1 * (AUT + Lincoln + Otago + Waikato)^-2 * (AUT +
## Lincoln + Waikato)^-4 * (AUT + Massey)^-4 * (AUT + Massey + Otago + Waikato)^-1 * (AUT +
## Massey + Victoria)^-1 * (AUT + Massey + Victoria + Waikato)^-3 * (AUT + Massey +
## Waikato)^-4 * (AUT + Otago)^-1 * (AUT + Otago + Waikato)^-1 * (AUT + Waikato)^-5 *
## Canterbury^40 * (Canterbury + Lincoln + Massey + Otago + UoA + Victoria + Waikato)^-7 *
## (Canterbury + Lincoln + Massey + UoA + Victoria + Waikato)^-2 * (Canterbury + Lincoln +
## Massey + Victoria + Waikato)^-2 * (Canterbury + Lincoln + Otago + UoA + Victoria +
## Waikato)^-2 * (Canterbury + Lincoln + Otago + Victoria + Waikato)^-2 * (Canterbury +
## Massey + Otago + UoA + Victoria + Waikato)^-2 * (Canterbury + Massey + Otago + Victoria +
## Waikato)^-1 * (Canterbury + Otago + UoA + Victoria + Waikato)^-1 * Lincoln^16 * (Lincoln
## + Massey)^-1 * (Lincoln + Massey + Otago + UoA + Victoria + Waikato)^-1 * (Lincoln +
## Massey + Otago + Victoria + Waikato)^-1 * (Lincoln + Massey + Victoria)^-2 * (Lincoln +
## Massey + Victoria + Waikato)^-3 * (Lincoln + Massey + Waikato)^-1 * (Lincoln + Otago +
## Victoria)^-1 * (Lincoln + Otago + Victoria + Waikato)^-2 * (Lincoln + Otago + Waikato)^-1
## * (Lincoln + Victoria)^-4 * (Lincoln + Victoria + Waikato)^-1 * (Lincoln + Waikato)^-1 *
## Massey^29 * (Massey + Otago + Victoria + Waikato)^-1 * Otago^40 * (Otago + UoA + Victoria
## + Waikato)^-1 * (Otago + Victoria + Waikato)^-2 * UoA^40 * Victoria^39 * (Victoria +
## Waikato)^-2 * Waikato^38)
Object universities is simply the support function for the
universities’ strengths for all years’ observations (assumed to be
independent) 2013-2020:
universities_maxp <- maxp(universities)
universities_maxp
## UoA Otago AUT Canterbury Lincoln Massey Victoria Waikato
## 0.6305043 0.1215205 0.0077547 0.0851247 0.0108839 0.0234650 0.0820183 0.0387286
pie(universities_maxp)
Note the dominance of Auckland (UoA). We observe that AUT has increased its strength since 2016 and we can test whether this increase is significant:
H1 <- hyper2(pnames=c("AUT_pre",universities_table[c(1:2,4:8),1]))
for(i in 2012:2016){
b <- f(universities_table,i)
rownames(b)[3] <- "AUT_pre"
H1 <- H1 + ordertable2supp(b)
} # column loop closes
H2 <- hyper2(pnames=c("AUT_post",universities_table[c(1:2,4:8),1]))
for(i in 2017:2020){
jj <- f(universities_table,i)
rownames(jj)[3] <- "AUT_post"
H2 <- H2 + ordertable2supp(jj)
}
maxp(H2)
## AUT_post UoA Otago Canterbury Lincoln Massey Victoria Waikato
## 0.015487 0.534311 0.139502 0.099503 0.046538 0.045748 0.075663 0.043247
H <- H1 + H2
maxp(H)
## AUT_post AUT_pre Canterbury Lincoln Massey Otago UoA Victoria Waikato
## 0.0086605 0.0039696 0.0736179 0.0081125 0.0215735 0.1190954 0.6543450 0.0762467 0.0343788
So from the above we see that the estimated post-2016 strength of AUT is larger than the pre-2016 strength. Is this change significant?
samep.test(H,c("AUT_pre","AUT_post"))
##
## Constrained support maximization
##
## data: H
## null hypothesis: AUT_pre = AUT_post
## null estimate:
## AUT_post AUT_pre Canterbury Lincoln Massey Otago UoA Victoria Waikato
## 0.0059933 0.0059933 0.0736396 0.0082884 0.0213885 0.1190810 0.6552033 0.0761355 0.0342771
## (argmax, constrained optimization)
## Support for null: -326.15 + K
##
## alternative hypothesis: sum p_i=1
## alternative estimate:
## AUT_post AUT_pre Canterbury Lincoln Massey Otago UoA Victoria Waikato
## 0.0086605 0.0039696 0.0736179 0.0081125 0.0215735 0.1190954 0.6543450 0.0762467 0.0343788
## (argmax, free optimization)
## Support for alternative: -325.34 + K
##
## degrees of freedom: 1
## support difference = 0.80926
## p-value: 0.2033
Following lines create universities.rda, residing in the data/
directory of the package.
save(universities_table,universities_maxp,universities,file="universities.rda")
hyper2 Package: Likelihood Functions for Generalized Bradley-Terry Models.” The R Journal 9 (2): 429–39.