Cuzick test

読んだ。
N Engl J Med. 2014 Jun 18.
N Engl J Med. 2014 Jun 18.
どちらの論文も、APOC3というアポリポ蛋白の話で、これの機能欠失型変異ではTGが低めでHDLが高め、結果として血管病変や心疾患が減る…
その中の解析でCuzike testというものを知った(Stat Med. 1985 Jan-Mar;4(1):87-90.)。
Wilcoxon検定の拡張で、群に傾向性があるかないかを検定する。ヨンクヒール・タピストラ検定みたいな感じ。
ぐぐったら実装を見つけた。

cuzick <-  function(x,z,test.type=c("two.sided", "upper", "lower")) {
	N = length(z)
	n = unique(z)
	ranks=rank(x)
	T = sum(ranks*z)
	p = (table(z)/N)
	E_Z = sum(unique(z)*p)
	E_T = 0.5*N*(N+1)*E_Z
	Var_Z = sum(unique(z)^2*p) - E_Z^2
	Var_T = N^2*(N+1)/12*Var_Z
	Zscore = (T-E_T)/sqrt(Var_T)
	if(test.type == "two.sided") {
		pval = 2*pnorm(-abs(Zscore))
	} else if(test.type == "upper") {
		pval = pnorm(Zscore,lower.tail=F)
	} else 	pval = pnorm(Zscore,lower.tail=T)
	out = data.frame(cbind(Zscore,pval,test.type))
	colnames(out) = c("Z","p","testType")
	return(out)
}

z <-  unlist(mapply(rep, 1:5, c(8, 10, 9, 9, 9)))
x <-  c(0, 0, 1, 1, 2, 2, 4, 9,
		0, 0, 5, 7, 8, 11, 13, 23, 25, 97,
		2, 3, 6, 9, 10, 11, 11, 12, 21,
		0, 3, 5, 6, 10, 19, 56, 100, 132,
		2, 4, 6, 6, 6, 7, 18, 39, 60)
cuzick(x, z, "two.sided")
boxplot(x ~ z)

                 Z                 p  testType
1 2.11038574142394 0.034825143619734 two.sided

微妙にありそうな…