メタアナリシス: 要約オッズ比を計算してみよう．

# Program for meta-analysis: Summary Odds Ratio.
# Awk script written by M.Yanagi.
# Only fixed-effects model.
#
#	Outcome		+	-	total
#	Factor	+	a	b	g
#		-	c	d	h
#	Total		e	f	n
#
# Data should be as	"a	b	c	d".
#
BEGIN	{
	FS = "\t"		# Separators are tabs.
	split(ARGV[1], k, ".")
	outfile = k[1]".out"
}
NR == 1 {
	print $0 > outfile	# The first line should be a comment.
	print "--- Odds Ratio and 95% CI for Each Trial ---" > outfile
	print "No.\ta\tb\tc\td\tOR\tMH_l_i\tMH_u_i\tX2" > outfile
}
NR>1	{
	i++
	printf("Reading  %d,\r", NR)
	a = $1;  b = $2;  c = $3;  d = $4;
	e = a+c; f = b+d; g = a+b; h = c+d; n = g+h;
	if (e*f*g*h == 0) { next }
	else if (a*b*c*d == 0) {
		a+=0.5; b+=0.5; c+=0.5; d+=0.5
		e+=1; f+=1; g+=1; h+=1; n+=2
	}
#
	OR         = a*d/(b*c)
	ORi[i]     = OR
	X2         = (a*d-b*c)**2*n/(e*f*g*h)
# ****** Mantel-Haenszel method ******
	Var_mh     = n/(b*c)
	Sum_w_OR  += OR/Var_mh
	Sum_w     += 1/Var_mh
	ORmh       = Sum_w_OR/Sum_w
#
	SumF      += F = a*d*(a+d)/ n**2
	SumG      += G = (a*d*(b+c) +  b*c*(a+d))/ n**2
	SumH      += H = b*c*(b+c)/ n**2
	SumR      += R = a*d/n
	SumS      += S = b*c/n
	Var_mh     = SumF/(SumR**2*2)+SumG/(2*SumR*SumS)+SumH/(SumS**2*2)
	rVmh       = sqrt(Var_mh)
#
	rVmh_i     = se_ln_ORmh = sqrt(1/a + 1/b + 1/c + 1/d)
	weight[i]  = 1/(1/a + 1/b + 1/c + 1/d)		# by RevMan help
	ORmh_l_i   = OR*exp(-1.96*rVmh_i)	# lower bound
	ORmh_u_i   = OR*exp(+1.96*rVmh_i)	# upper bound
	ORmh_l     = ORmh*exp(-1.96*rVmh)
	ORmh_u     = ORmh*exp(+1.96*rVmh)
# ****** Peto method ******
	E          = e*g/n
	O_E        = a-E
	Var_p      = E*f*h/(n*(n-1))
	X2p        = O_E**2/Var_p
	Sum_O_E   += O_E
	Sum_Var_p += Var_p
	lnORp      = Sum_O_E/Sum_Var_p
	ORp        = exp(lnORp)
	ORp_l      = exp(lnORp - 1.96/sqrt(Sum_Var_p))
	ORp_u      = exp(lnORp + 1.96/sqrt(Sum_Var_p))
	Sum_X2p   += X2p
	Qp         = Sum_X2p-Sum_O_E**2/Sum_Var_p
						# Peto homogeneity
#
	printf("%d\t%d\t%d\t%d\t%d\t%7.4f %7.4f %7.4f %7.4f\n",
		i,a,b,c,d,OR,ORmh_l_i,ORmh_u_i,X2) > outfile
}
END	{
	print ""
	print "" > outfile
	print "------ Summary Odds Ratio ------" > outfile
	print "OR_MH\tMH_low\tMH_upp\tVarMH\t\tOR_Peto\tPeto_l\tPeto_u\tX2p"\
		> outfile
	printf("%7.4f %7.4f %7.4f %7.4f \t%7.4f %7.4f %7.4f %7.4f\n",
	  ORmh,ORmh_l,ORmh_u,Var_mh,ORp,ORp_l,ORp_u,Sum_X2p) > outfile
	print "------ Test of homogeneity ------" > outfile
	for(j=1;j<=i;j++) {
		Qmh += weight[j]*(( log(ORmh)-log(ORi[j]) )**2)
	}
	print "Peto Q\tMH Q\tfreedom = "i-1 > outfile
	printf("%7.4f %7.4f\n", Qp,Qmh) > outfile
	print "" > outfile
	print "X2@df 0.05 = (5.991, 7.815, 9.488, 11.07, 12.59, 14.0, 15.51, 16.92)@(2,3,4,5,6,7,8,9)" > outfile
}
# < Reference >
# 1) Petitti DB:  Meta-analysis, decision analysis, and cost-effective
#    analysis.  Oxford, Oxford University Press, 1994.
# 2) Hamajima N.  Journal of the Japanese Association for
#    Cerebro-cardiovascular Disease Control  33(3):266-271, 1998.
# 3) Cochrane collaboration:  Statistical methods programmed in metaview.
#    In: Review manager help version 4.1, 2000.
# "Where zero cells cause problems with computation of effects or standard
#  errors, 0.5 is added to all cells (ai, bi, ci, di) for that study, except
#  when ai = ci = 0 or bi = di = 0, when the relative effect measures ORi and
#  RRi are undefined."

Test data
6	22	19	142
5	17	10	68

Test data
--- Odds Ratio and 95% CI for Each Trial ---
No.	a	b	c	d	OR	MH_l_i	MH_u_i	X2
1	6	22	19	142	 2.0383  0.7337  5.6629  1.9261
2	5	17	10	68	 2.0000  0.6037  6.6263  1.3209

------ Summary Odds Ratio ------
OR_MH	MH_low	MH_upp	VarMH	OR_Peto	Peto_l	Peto_u	X2p
 2.0216  0.9289  4.4001  0.1574 	 2.2375  0.9282  5.3936  3.2236
------ Test of homogeneity ------
Peto Q	MH Q	freedom = 1
 0.0052  0.0006

X2@df 0.05 = (5.991, 7.815, 9.488, 11.07, 12.59, 14.0, 15.51, 16.92)@(2,3,4,5,6,7,8,9)

• 各試験毎のオッズ比とその95％信頼区間，カイ2乗値．
• 要約オッズ比の値とその95％信頼区間（Mantel-Haenszel法とPeto法）．
• 均一性テストの結果．

リスク比の要約（NNTの計算のために）