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

廿TT

譬如水怙牛過窓櫺 頭角四蹄都過了 因甚麼尾巴過不得

指数分布に従う確率変数の和の分布をたたみこみで求める

確率分布 R

計算

X をパラメータ λ の指数分布に従う確率変数, Y をパラメータ μ の指数分布に従う確率変数とする.

それぞれの密度関数を f_X(x), f_Y(y) と表記する.

XY は互いに独立とする.

Z = X + Y と置くと XY (= Z - X) の結合密度関数は

\displaystyle f(x,z-x) = \lambda e^{-\lambda x} \mu e^{-\mu (z-x)}\\
\displaystyle = \lambda \mu e^{-\mu z} e^{(\mu -\lambda)x}.

Z は以下の密度関数を持つ.

\displaystyle f_Z(z) = \int^{z}_{0} f_X(x) \, f_Y(z-x) \,dx \\
\displaystyle = \lambda \mu e^{-\mu z}\int^{z}_{0}e^{(\mu -\lambda)x} \,dx \\
\displaystyle = \lambda \mu e^{-\mu z}\cdot \frac{1}{\mu - \lambda} \left[ e^{(\mu -\lambda)x}\right]^{z}_{0} \\
\displaystyle = \lambda \mu e^{-\mu z}\cdot \frac{1}{\mu - \lambda}\left\{e^{(\mu -\lambda)z}-1 \right\} \\
\displaystyle = \frac{\lambda \mu}{\mu - \lambda}(e^{-\lambda z} - e^{-\mu z}).

分布関数は,

\displaystyle F_Z(z) = \int^{z}_{0} f_Z(t) \, dt \\
\displaystyle = \frac{\lambda \mu}{\mu - \lambda}\int^{z}_{0} (e^{-\lambda t} - e^{-\mu t})\, dt \\
\displaystyle = \frac{\lambda \mu}{\mu - \lambda} \left[- \frac{1}{\lambda}e^{-\lambda t} + \frac{1}{\mu}e^{-\mu t} \right]^{z}_{0} \\
\displaystyle =\frac{\lambda \mu}{\mu - \lambda}\left\{ -\frac{1}{\lambda}(e^{-\lambda z}-1) + \frac{1}{\mu}(e^{-\mu z} -1) \right\}  \\
\displaystyle = - \frac{\mu}{\mu -\lambda} (e^{-\lambda z} -1) + \frac{\lambda}{\mu -\lambda} (e^{-\mu z} -1) \\
\displaystyle = 1 - \frac{\mu}{\mu -\lambda} e^{- \lambda z}+ \frac{\lambda}{\mu -\lambda} e^{- \mu z}.

シミュレーション

密度関数

パラメータの異なる指数乱数をふたつ発生させ, 足し合わせて分布を見る.

上で求めた密度関数をヒストグラムに重ねてプロットした.

#R のコード
lambda <- 2
mu <- 8
z <- rexp(10000,lambda) + rexp(10000,mu)
hist(z,freq=FALSE,breaks="scott")
dsexp <- function(z){
  ((lambda*mu)/(mu-lambda))*(exp(-lambda*z)-exp(-mu*z))
}
curve(dsexp(x),lwd=1.5,col="red3",add=TRUE)

f:id:abrahamcow:20150523195642p:plain

分布関数

同様, 経験分布関数に上述の分布関数を重ねる.

psexp <- function(z){
1 - (mu/(mu-lambda))*exp(-lambda*z) +
(lambda/(mu-lambda))*exp(-mu*z)
}
plot.ecdf(z, verticals= TRUE, do.p = FALSE,main="",xlab="z",ylab="Fn(z)")
curve(psexp(x),add=TRUE,col="red3",from=0)

f:id:abrahamcow:20150523195736p:plain

リパラメタライゼーション

X をパラメータ 1/λ の指数分布に従う確率変数, Y をパラメータ 1/μ の指数分布に従う確率変数とすると,

\displaystyle f_Z(z) = \frac{1}{\lambda-\mu}(e^{-\frac{1}{\lambda} z} - e^{-\frac{1}{\mu} z}).

\displaystyle F_Z(z) =  1 - \frac{\lambda}{\lambda-\mu}e^{-\frac{1}{\lambda}z} +\frac{\mu}{\lambda-\mu}e^{-\frac{1}{\mu}z}

lambda <- 8
mu <- 2
z <- rweibull(10000,1,lambda) + rweibull(10000,1,mu)
#z <- rexp(10000,1/lambda) + rexp(10000,1/mu)
hist(z,breaks="scott",freq=FALSE)
dsexp <- function(z){
  (exp(-z/lambda)-exp(-z/mu))/(lambda-mu)
}
curve(dsexp(x),col="red3",add=TRUE)

f:id:abrahamcow:20150626233420p:plain

psexp <- function(z){
  1 - (lambda/(lambda-mu))*exp(-z/lambda) +
    (mu/(lambda-mu))*exp(-z/mu)
}
plot.ecdf(z, verticals= TRUE, do.p = FALSE,main="",xlab="z",ylab="Fn(z)")
curve(psexp(x),add=TRUE,col="red3",from=0)

f:id:abrahamcow:20150626233425p:plain

関連エントリ

abrahamcow.hatenablog.com
abrahamcow.hatenablog.com