読者です 読者をやめる 読者になる 読者になる

INFINITY -数学とかプログラミングとか-

統計とプログラムを使って役に立たせたい

TeX用コマンド入力を支援するための辞書をご利用ください。
sanctuary's blogは,適当なことが書いてあります。

双方な分布のネゲントロピーに関する考察

分布3\mu \phi(3x) + 3(1 - \mu)\phi(3(x - 1))を平均0,分散1にしたときのネゲントロピーを計算し,横軸が0.001刻みで\muを0から1まで動かし,縦軸がネゲントロピーの近似した図
(ちゃんとやったのに図が間違っているかも?)

f:id:infinity_th4:20121022175226p:plain


ネゲントロピーが大きいという事はガウス分布から形状が離れているという事である。そこでどういう意味で形状が離れているかという話で,最も分布が双方になった時,つまり,\mu=0.5のときに最大であって欲しかったわけだが,そうではない(らしい)。(図が間違っているのかも?

\mu=0.2 or 0.8あたりで最大なので,対称性からのズレか?

間違っていると思うので,誰か教えてください(^_^;)

  • コード

#mu*3*phi(3*x) + (1 - mu)*3*phi(3*(x - 1))を平均0分散1にしたときのネゲントロピーの計算
phi=function(x) dnorm(x,mean=0,sd=1)
p = function(x){
    mu*3*phi(3*x) + (1 - mu)*3*phi(3*(x - 1))
    
}

mu=0
mean=1-mu
var=((10-9*mu)/9)-(1-mu)^2

pdf=function(x){
    p(sqrt(var)*x+mean)*(sqrt(var))
}

integrate(pdf,-Inf,Inf)

#双方な分布のグラフ
xi=seq(-50,50,length=1000)

for( mu in 0:10){
    mu = mu*0.1
    mean=1-mu
    var=((10-9*mu)/9)-(1-mu)^2
    plot(xi,pdf(xi),xlab="xi",xlim=c(-5,3),ylim=c(0,0.8),type="l",col=mu*10,lwd=1)
    par(new=T)
}
par(new=T)
plot(xi,phi(xi),xlab="xi",xlim=c(-5,3),ylim=c(0,0.8),type="o",col=1,lwd=2)
par(new=T)

mu = 0.2
mean=1-mu
var=((10-9*mu)/9)-(1-mu)^2
plot(xi,pdf(xi),xlab="xi",xlim=c(-5,3),ylim=c(0,0.8),type="l",col=mu*10,lwd=2)
par(new=T)
mu = 0.8
mean=1-mu
var=((10-9*mu)/9)-(1-mu)^2
plot(xi,pdf(xi),xlab="xi",xlim=c(-5,3),ylim=c(0,0.8),type="l",col=mu*10,lwd=3)
par(new=T)
mu = 0.5
mean=1-mu
var=((10-9*mu)/9)-(1-mu)^2
plot(xi,pdf(xi),xlab="xi",xlim=c(-5,3),ylim=c(0,0.8),type="l",col=mu*10,lwd=4)






e1=function(x){
    x*pdf(x)
}
e2=function(x){
    x^2*pdf(x)
}
integrate(e1,-Inf,Inf)
integrate(e2,-Inf,Inf)


#gauss=function(x) dnorm(x,mean=1-mu,sd=sqrt((1+7*mu-4*mu^2)/4))



k1=36/(8*sqrt(3)-9)
k2a=1/(2-6/pi)

#g1:奇関数
g1=function(x){
    x*exp(-(x^2)/2)*pdf(x)
}
#g2a:偶関数
g2a=function(x){
    abs(x)*pdf(x)
}

integrate(g1,-Inf,Inf)$val

integrate(g2a,-Inf,Inf)$val




mu=0

mean=1-mu
var=((10-9*mu)/9)-(1-mu)^2
Ja(mu)

mu=1
mean=1-mu
var=((10-9*mu)/9)-(1-mu)^2
Ja(mu)

for(mu in 0:1000){
    mu = mu*0.001
    mean=1-mu
    var=((10-9*mu)/9)-(1-mu)^2
    plot(mu,Ja(mu),xlab="mu",xlim=c(0,1),ylim=c(0,0.08),type="p",lwd=1)
    par(new=T)
    
}