# 2005 年 1 月 18 日作成
# Large Wishart matrix の固有根の経験分布と極限分布(Winger law)の同時の作図
#
#
# 倍数
ratio<-1.2
y<-1/ratio
# 次元
no.of.dim<-100
# 標本数
no.of.sample<-trunc(ratio*no.of.dim)
print(no.of.sample)
emp<-rep(0,no.of.dim)
wishart.matrix<-matrix(rep(0,no.of.dim**2),no.of.dim)
#print(wishart.matrix)
for (i in 1:no.of.sample){
x<-rnorm(no.of.dim,0,1)
wishart.matrix<-wishart.matrix+x%*%t(x)/no.of.sample
#print(wishart.matrix)
}
# Wishart matrix の固有根
emp<-c(eigen(wishart.matrix)$values)
# 理論的な下端と上端
a<-1+y-2*sqrt(y)
b<-1+y+2*sqrt(y)
# 固有根の経験分布の作図
plot(sort(emp),(1:no.of.dim)/no.of.dim,type="s",xlim=c(0,b+1),ylim=c(0,1),main="emp and Winger (p=100)")
par(new=T)
winger<-function(x){
return(sqrt((x-a)*(b-x))/(2*pi*y*x))
}
mesh<-0.1
k<-trunc((b-a)/mesh)+1
x<-seq(a,b,by=mesh)
cwinger<-rep(0,k)
for(i in 1:k){
cwinger[i]<-integrate(winger,a,x[i])
}
plot(x,cwinger,xlim=c(0,b+1),ylim=c(0,1),type="l",xlab="",ylab="")