(4) ケース・コントロール比が大きいとき(コホート疫学)

# コホート
# 2x2分割表
NOp<-FALSE
N1<-10000
N2<-20

a<-1000
A<-c(rep(1,a),rep(2,N1-a))
A<-sample(A)
xlim<-ylim<-c(0,1)

bs<-seq(from=1,to=N2-1,by=1)
t<-seq(from=0,to=1,length=100)*2*pi
ra<-0.2
rb<-ra
max<-0.25
may<-0.5
mbx<-0.75
mby<-0.5

p.out<-rep(0,length(bs))

for(i in 1:length(bs)){
	b<-bs[i]
	B<-c(rep(1,b),rep(2,N2-b))
	B<-sample(B)
	plot(ra*cos(t)+max,ra*sin(t)+may,type="l",xlim=xlim,ylim=ylim,main="",xlab="",ylab="")
	par(new=TRUE)
	plot(rb*cos(t)+mbx,rb*sin(t)+mby,type="l",xlim=xlim,ylim=ylim,main="",xlab="",ylab="")
	par(new=TRUE)
	rs<-runif(N1)*ra*0.9
	ts<-runif(N1)*2*pi
	plot(rs*cos(ts)+max,rs*sin(ts)+may,pch=19,cex=0.1,col=A,xlim=xlim,ylim=ylim,main="",xlab="",ylab="")
	par(new=TRUE)
	rs<-runif(N2)*rb*0.9
	ts<-runif(N2)*2*pi
	#print(table(as.data.frame(cbind(A,B))))
	v1<-c(rep(1,length(A)),rep(2,length(B)))
	v2<-c(A,B)
	print(table(v1,v2))
	p.f<-fisher.test(table(v1,v2))$p.value
	p.out[i]<-p.f
	if(NOp)p.f=""
	plot(rs*cos(ts)+mbx,rs*sin(ts)+mby,pch=19,cex=1,col=B,xlim=xlim,ylim=ylim,main=p.f,xlab="",ylab="")


	par(ask=TRUE)
}

plot(bs,p.out)