🦢

【ワイブル分布】補足:打ち切りする場合の形状パラメータの最尤推定量の期待値(総数 n → ∞)

2024/08/17に公開

はじめに

今回はワイブル分布の形状パラメータ m の最尤推定量 \hat{m} を不偏推定量(期待値がパラメータの真値と等しい推定量)とする補正係数に関する補足です。全体で n 個のサンプルについて試験を同時に実施し、 r 個の観測値が得られた段階で打ち切る場合を考えます。補正係数を求めるために、 最尤推定量 \hat{m} の期待値を求める必要がありますが、解析的に求めることができません。本記事では、補正係数の近似式の立式のため、 r を一定として、 n を無限大とした場合の最尤推定量 \hat{m} の期待値が \dfrac{r}{r-2}m となることを導出します。導出には変数変換を多用しますので、統計検定の勉強において確率変数の変換の復習となるかもしれません。
形状パラメータ m の補正係数の詳細についてはこちら
ワイブル分布の最尤推定量に関する詳細についてはこちら
補正係数を用いた最尤推定量の改良方法についてはこちらをご覧ください。

概要

以下、前置きです。

  • 確率変数 X_{1},\ ...\ ,X_{n} は互いに独立に形状パラメータ m 、尺度パラメータ \eta のワイブル分布 W(m,\eta) に従い、その確率密度関数を f(x)、累積分布関数を F(x) とする
  • X_{(1)},\ ...\ ,X_{(n)}X_{1},\ ...\ ,X_{n} の順序統計量とする
    \big(X_{(i)}X_{1},\ ...\ ,X_{n} の中で i 番目に小さい確率変数\big)
  • X_{i} から得られる観測値を x_{i}X_{(i)} から得られる観測値を x_{(i)} とする(i=1,\ ...\ ,n)
  • 打切りについては n 個のサンプルの試験において、小さい順に r 個の観測値 x_{(1)},\ ...\ ,x_{(r)} が得られた段階で試験を打ち切ると仮定

また、確率密度関数 f(x) 、累積分布関数 F(x)

\begin{aligned} f(x)&=\dfrac{m}{\eta}\left(\dfrac{x}{\eta} \right)^{m-1} \exp\left\{-\left(\dfrac{x}{\eta} \right)^{m} \right\} \\ F(x)&=1-\exp\left\{-\left(\dfrac{x}{\eta} \right)^{m} \right\} \\ \end{aligned}

と表されます。

小さい順に r 個の観測値 x_{(1)},\ ...\ ,x_{(r)} が得られた段階で試験を打ち切るとした場合の形状パラメータ m の最尤推定量 \hat{m}_{X} は以下の式を満たします。

\begin{aligned} \dfrac{1}{\hat{m}_{X}} + \dfrac{1}{r}\sum_{i=1}^{r}\ln{x_{(i)}} -\dfrac{\sum_{i=1}^{r} x_{(i)}^{\hat{m}_{X}} \ln{x_{(i)}} +(n-r) x_{(r)}^{\hat{m}_{X}}\ln x_{(r)}} {\sum_{i=1}^{r} x_{(i)}^{\hat{m}_{X}} +(n-r) x_{(r)}^{\hat{m}_{X}} } =0 \end{aligned}

このとき、r を一定として n を無限大としたときの最尤推定量の期待値 E[\hat{m}_{X}] は以下の式で表すことができます。

\begin{aligned} \lim_{n\to \infty}E[\hat{m}_{X}]=\dfrac{r}{r-2}m \end{aligned}

k_{r,n}=\dfrac{\hat{m_{X}}}{m} とすると、以下のように表せます。

\begin{aligned} \lim_{n\to \infty}E[k_{r,n}]=\dfrac{r}{r-2} \end{aligned}

導出については以下の手順で行います。

  1. r 個の観測値が得られた段階で試験を打ち切るとした場合の形状パラメータ m の最尤推定量 を求める
  2. n を無限大としたときの最尤推定量の期待値を求める 
    期待値は積分の形で表し、以下の順で変数変換
    2.1. Z_{(i)} = \left( \dfrac{X_{(i)}}{\eta} \right)^{m} として X_{(i)} \to Z_{(i)} の変数変換
        \to Z_{(i)} \sim Exp(1) の順序統計量
     
    2.2. W_{i}= (n-i+1)(Z_{(i)}-Z_{(i-1)}) として Z_{(i)} \to W_{i} の変数変換 (Z_{(0)}=0)
        \to W_{i} \sim Exp(1) (互いに独立)
     
    2.3. \displaystyle \lim_{n \to \infty}\hat m_{X} を求め、E[\hat{m}_{X}] \to\displaystyle \lim_{n \to \infty}E[\hat m_{X}] の形に変形し、n を除去
     
    2.4. V=W_{1}+ \cdots +W_{r},\ U_{(i)} = \dfrac{\sum_{j=1}^{i} W_{j}}{V}\hspace{10pt} (i=1,\ ...\ ,r-1) として、
       W_{1},\ ...\ ,W_{r} \to U_{(1)},\ ... \ ,U_{(r-1)}, V の変数変換
       \to U_{(i)} (1,\ ...\ ,r-1) \sim 標準一様分布 U(0,1) の順序統計量、 V \sim Gamma(r,1)
     
    2.5. T_{i} = -\ln U_{(i)}\ (i=2,\ ...\ ,r-1),\ S= - \ln(U_{(1)} \cdots U_{(r-1)}) として、
       U_{(1)},\ ...\ U_{(r-1)} \to T_{2},\ ...\ ,T_{r-1},S の変数変換

1. 打ち切りする場合の形状パラメータ m の最尤推定量

まず、 n 個の部品の寿命試験において、小さい順に r 個の観測値 x_{(1)},\ ...\ ,x_{(r)} が得られた段階で試験を打ち切る場合の形状パラメータ m の最尤推定量を求めます。残りの n-r 個については観測値は得られませんが、小さい順に r 番目の観測値 x_{(r)}より大きいという条件があります。

この場合の尤度関数 L(m,\eta)

\begin{aligned} L(m,\eta) &= P(X_{(1)}=x_{(1)},\ ...\ ,X_{(r)}=x_{(r)},\ X_{(r+1)},\ ...\ ,X_{(n)} \gt x_{(r)}) \\ &={}_{n}P_{r} \left(\prod_{i=1}^{r} f(x_{(i)}) \right) \{1-F(x_{(r)})\}^{n-r} \\ &={}_{n}P_{r} \left(\prod_{i=1}^{r} \dfrac{m}{\eta^{m}}x_{(i)}^{m-1} \exp\left\{-\left( \dfrac{x_{(i)}}{\eta} \right)^{m}\right\} \right) \exp\left\{-(n-r)\left( \dfrac{x_{(r)}}{\eta} \right)^{m}\right\} \\ &={}_{n}P_{r} \dfrac{m^{r}}{\eta^{mr}} \left(\prod_{i=1}^{r} x_{(i)}^{m-1} \right) \exp\left\{-\dfrac{1}{\eta^{m}}\left(\sum_{i=1}^{r}x_{(i)}^{m} + (n-r)x_{(r)}^{m} \right)\right\} \\ \end{aligned}

と表せます。対数尤度関数 l(m,\eta) については

\begin{aligned} l(m,\eta) &=\ln L(m,\eta) \\ &=\ln {}_{n}P_{r} +r\ln m -mr\ln \eta +(m-1)\sum_{i=1}^{r} \ln x_{(i)} -\dfrac{1}{\eta^{m}}\left(\sum_{i=1}^{r}x_{(i)}^{m} + (n-r)x_{(r)}^{m} \right) \\ \end{aligned}

対数尤度関数が最大となる m,\eta を求めると、以下の2式が成り立ちます。

\left\{\begin{aligned} \dfrac{\partial}{\partial m}\ln{L(m,\eta)} = 0 \\ \dfrac{\partial}{\partial \eta}\ln{L(m,\eta)} = 0 \\ \end{aligned}\right.

\eta については

\begin{aligned} \dfrac{\partial}{\partial \eta}\ln{L(m,\eta)} &=-\dfrac{mr}{\eta} +\dfrac{m}{\eta^{m+1}} \left(\sum_{i=1}^{r}x_{(i)}^{m} + (n-r)x_{(r)}^{m} \right) \\ &=\dfrac{m}{\eta} \left(-r + \dfrac{\sum_{i=1}^{r}x_{(i)}^{m} + (n-r)x_{(r)}^{m}}{\eta^{m}} \right) =0\\ \end{aligned}

よって、

\begin{aligned} r &=\dfrac{\sum_{i=1}^{r}x_{(i)}^{m} + (n-r)x_{(r)}^{m}}{\eta^{m}} \\ \eta^{m} &= \dfrac{1}{r} \left( \sum_{i=1}^{r}x_{(i)}^{m} + (n-r)x_{(r)}^{m} \right) \\ \eta &= \sqrt[\footnotesize m]{\dfrac{1}{r} \left( \sum_{i=1}^{r}x_{(i)}^{m} + (n-r)x_{(r)}^{m} \right)} \end{aligned}

m については

\begin{aligned} \dfrac{\partial}{\partial m}\ln{L(m,\eta)} &=\dfrac{r}{m} -r\ln \eta + \sum_{i=1}^{r}\ln{x_{(i)}} -\left\{ \sum_{i=1}^{r} \left( \dfrac{x_{(i)}}{\eta}\right)^{m} \ln \dfrac{x_{(i)}}{\eta} + (n-r)\left(\dfrac{x_{(r)}}{\eta} \right)^{m}\ln \dfrac{x_{(r)}}{\eta} \right\}\\ &=\dfrac{r}{m} -r\ln \eta + \sum_{i=1}^{r}\ln{x_{(i)}} -\dfrac{1}{\eta^{m}}\left\{ \sum_{i=1}^{r} x_{(i)}^{m} \ln x_{(i)} + (n-r)x_{(r)}^{m} \ln x_{(r)} -\sum_{i=1}^{r} x_{(i)}^{m} \ln \eta - (n-r)x_{(r)}^{m} \ln \eta \right\} \\ &=\dfrac{r}{m} -r\ln \eta + \dfrac{\sum_{i=1}^{r} x_{(i)}^{m} + (n-r)x_{(r)}^{m} }{\eta^{m}} \ln \eta + \sum_{i=1}^{r}\ln{x_{(i)}} -\dfrac{1}{\eta^{m}}\left\{ \sum_{i=1}^{r} x_{(i)}^{m} \ln x_{(i)} + (n-r)x_{(r)}^{m} \ln x_{(r)} \right\}\\ &=\dfrac{r}{m} + \sum_{i=1}^{r}\ln{x_{(i)}} -r\dfrac{\sum_{i=1}^{r} x_{(i)}^{m} \ln x_{(i)} + (n-r)x_{(r)}^{m} \ln x_{(r)}} {\sum_{i=1}^{r}x_{(i)}^{m} + (n-r)x_{(r)}^{m}} \\ &=r\left( \dfrac{1}{m} + \dfrac{1}{r}\sum_{i=1}^{r}\ln{x_{(i)}} -\dfrac{\sum_{i=1}^{r} x_{(i)}^{m} \ln x_{(i)} + (n-r)x_{(r)}^{m} \ln x_{(r)}} {\sum_{i=1}^{r}x_{(i)}^{m} + (n-r)x_{(r)}^{m}} \right) = 0\\ \end{aligned}

よって、確率変数 X_{(i)} についての最尤推定量 \hat m_{X} は以下の式を満たします。

\begin{aligned} \dfrac{1}{\hat m_{X}} + \dfrac{1}{r}\sum_{i=1}^{r}\ln{x_{(i)}} -\dfrac{\sum_{i=1}^{r} x_{(i)}^{\hat m_{X}} \ln x_{(i)} + (n-r)x_{(r)}^{\hat m_{X}} \ln x_{(r)}} {\sum_{i=1}^{r}x_{(i)}^{\hat m_{X}} + (n-r)x_{(r)}^{\hat m_{X}}} = 0 \\ \end{aligned}

となります。

2. n を無限大としたときの最尤推定量の期待値

次に、最尤推定量 \hat m_{X} の期待値 E[\hat m_{X}] を求めます。最尤推定量は解析的に解くことができないため、\hat m_{X} を使用して積分の形で表します。
導出には変数変換を多用します。導出までの流れを以下に記載します。
特記なき場合は i=1,\ ...\ ,r とします。

はじめに、X_{(i)} は形状パラメータ m 、 尺度パラメータ \eta のワイブル分布 W(m,\eta) の順序統計量とします。
 
2.1. Z_{(i)} = \left( \dfrac{X_{(i)}}{\eta} \right)^{m} として X_{(i)} \to Z_{(i)} の変数変換
    \to Z_{(i)} \sim Exp(1) の順序統計量
 
2.2. W_{i}= (n-i+1)(Z_{(i)}-Z_{(i-1)}) として Z_{(i)} \to W_{i} の変数変換 (Z_{(0)}=0)
    \to W_{i} \sim Exp(1) (互いに独立)
 
2.3. \displaystyle \lim_{n \to \infty}\hat m_{X} を求め、E[\hat{m}_{X}] \to\displaystyle \lim_{n \to \infty}E[\hat m_{X}] の形に変形し、n を除去
 
2.4. V=W_{1}+ \cdots +W_{r},\ U_{(i)} = \dfrac{\sum_{j=1}^{i} W_{j}}{V}\hspace{10pt} (i=1,\ ...\ ,r-1) として、
   W_{1},\ ...\ ,W_{r} \to U_{(1)},\ ... \ ,U_{(r-1)}, V の変数変換
   \to U_{(i)} (1,\ ...\ ,r-1) \sim 標準一様分布 U(0,1) の順序統計量、 V \sim Gamma(r,1)
 
2.5. T_{i} = -\ln U_{(i)}\ (i=2,\ ...\ ,r-1),\ S= - \ln(U_{(1)} \cdots U_{(r-1)}) として、
   U_{(1)},\ ...\ U_{(r-1)} \to T_{2},\ ...\ ,T_{r-1},S の変数変換

最尤推定量の \hat{m}_{X} 期待値 E[\hat{m}_{X}] については、

\begin{aligned} E[\hat{m}_{X}] &=\int_{0}^{\infty} \cdots \int_{0}^{x_{(2)}} \hat{m}_{X}\ L(m,\eta) \ dx_{(1)} \cdots dx_{(r)} \\ &=\int_{0}^{\infty} \cdots \int_{0}^{x_{(2)}} \hat{m}_{X}\ {}_{n}P_{r} \left(\prod_{i=1}^{r} \dfrac{mx_{(i)}^{m-1}}{\eta^{m}} \exp\left\{-\left( \dfrac{x_{(i)}}{\eta} \right)^{m}\right\} \right) \exp\left\{-(n-r)\left( \dfrac{x_{(r)}}{\eta} \right)^{m}\right\} dx_{(1)} \cdots dx_{(r)} \\ &=\int_{0}^{\infty} \cdots \int_{0}^{x_{(2)}} \hat{m}_{X}\ {}_{n}P_{r} \left(\prod_{i=1}^{r} \dfrac{mx_{(i)}^{m-1}}{\eta^{m}} \right) \exp\left\{-\dfrac{1}{\eta^{m}}\left(\sum_{i=1}^{r}x_{(i)}^{m} + (n-r)x_{(r)}^{m} \right)\right\} dx_{(1)} \cdots dx_{(r)} \\ \end{aligned}

と表せます。

2.1. X(i) から Z(i) の変数変換

まず、Z_{(i)}= \left( \dfrac{X_{(i)}}{\eta} \right)^{m} として X_{(i)} \to Z_{(i)} と変数変換します。

z_{(i)} = \left( \dfrac{x_{(i)}}{\eta} \right)^{m} (i=1,\ ...\ ,r) とすると、

\begin{aligned} J&=\begin{vmatrix} \dfrac{\partial x_{(1)}}{\partial z_{(1)}} & \cdots & \dfrac{\partial x_{(r)}}{\partial z_{(1)}} \\ \vdots & \ddots & \vdots \\ \dfrac{\partial x_{(1)}}{\partial z_{(r)}} & \cdots & \dfrac{\partial x_{(r)}}{\partial z_{(r)}} \\ \end{vmatrix} \\ &=\dfrac{\partial x_{(1)}}{\partial z_{(1)}} \cdots \dfrac{\partial x_{(r)}}{\partial z_{(r)}} \\ &=\dfrac{1}{\dfrac{\partial z_{(1)}}{\partial x_{(1)}} \cdots \dfrac{\partial z_{(r)}}{\partial x_{(r)}}} \\ &=\dfrac{1}{ \dfrac{mx_{(1)}^{m-1}}{\eta^{m}} \cdots \dfrac{mx_{(r)}^{m-1}}{\eta^{m}} } \\ &=\dfrac{\eta^{mr}}{ mx_{(1)}^{m-1} \cdots mx_{(r)}^{m-1}} \\ &=\left( \prod_{i=1}^{r} \dfrac{\eta^{m}}{ mx_{(i)}^{m-1} } \right)\\ \end{aligned}

となることから、期待値は以下のようになります。

\begin{aligned} E[\hat{m_{X}}] &=\int_{0}^{\infty} \cdots \int_{0}^{z_{(2)}} \hat{m}_{Z}\ {}_{n}P_{r} \exp\left\{-\left(\sum_{i=1}^{r}z_{(i)} + (n-r)z_{(r)} \right)\right\} dz_{(1)} \cdots dz_{(r)} \\ \end{aligned}

2.2. Z(i) から Wi の変数変換

次に W_{i} を以下のように定義し、Z_{(i)} \to W_{i} の変数変換をします。

\left\{\begin{aligned} W_{1} &= nZ_{(1)} \\ W_{2} &= (n-1) (Z_{(2)} - Z_{(1)}) \\ &\vdots \\ W_{r} &= (n-r+1) (Z_{(r)}- Z_{(r-1)}) \\ \end{aligned}\right.

とすると、

\left\{\begin{aligned} Z_{(1)} &= \dfrac{W_{1}}{n} \\ Z_{(2)} &= \dfrac{W_{1}}{n} + \dfrac{W_{2}}{n-1} \\ &\vdots \\ Z_{(r)} &= \sum_{i=1}^{r} \dfrac{W_{i}}{n-i+1} \\ \end{aligned}\right.
\begin{aligned} J&=\begin{vmatrix} \dfrac{\partial z_{(1)}}{\partial w_{1}} & \dfrac{\partial z_{(2)}}{\partial w_{1}} & \cdots & \dfrac{\partial z_{(r)}}{\partial w_{1}} \\ \dfrac{\partial z_{(1)}}{\partial w_{2}} & \dfrac{\partial z_{(2)}}{\partial w_{2}} & \cdots & \dfrac{\partial z_{(r)}}{\partial w_{2}} \\ \vdots & \vdots & \ddots & \vdots \\ \dfrac{\partial z_{(1)}}{\partial w_{r}} & \dfrac{\partial z_{(2)}}{\partial w_{r}} & \cdots & \dfrac{\partial z_{(r)}}{\partial w_{r}} \\ \end{vmatrix} \\ &=\begin{vmatrix} \dfrac{1}{n} & \dfrac{1}{n} & \cdots & \dfrac{1}{n} \\\\ 0 & \dfrac{1}{n-1} & \cdots & \dfrac{1}{n-1} \\\\ \vdots & \vdots & \ddots & \vdots \\\\ 0 & 0 & \cdots & \dfrac{1}{n-r+1} \\ \end{vmatrix} \\ &=\dfrac{1}{{}_{n}{P}_{r}} \\ \end{aligned}

また、

\begin{aligned} \sum_{i=1}^{r} W_{r} = \sum_{i=1}^{r} X_{(i)} +(n-r)X_{ (r)} \\ \end{aligned}

となることから E[\hat{m}_{X}]

\begin{aligned} E[\hat{m_{X}}] &=\int_{0}^{\infty} \cdots \int_{0}^{\infty} \hat{m}_{W} \exp\left( -\sum_{i=1}^{r}w_{i} \right) dw_{1} \cdots dw_{r} \\ \end{aligned}

と変形できます。

2.3. n を無限大とした場合の最尤推定量と期待値

2.2. による変数変換により、 n が存在するのは最尤推定量 \hat m_{W} のみになります。この段階で最尤推定量 \hat m_{W} について n \to \infty の極限をとり、 \displaystyle \lim_{n \to \infty} E[m_{X}] を求め、 n を除去します。

\displaystyle \lim_{n \to \infty} E[m_{X}] については

\begin{aligned} \lim_{n\to \infty }E[\hat{m_{X}}] &=\lim_{n\to \infty }\int_{0}^{\infty} \cdots \int_{0}^{\infty} \hat{m}_{W} \exp\left( -\sum_{i=1}^{r}w_{i} \right) dw_{1} \cdots dw_{r} \\ &=\int_{0}^{\infty} \cdots \int_{0}^{\infty} \lim_{n\to \infty } \hat{m}_{W} \exp\left( -\sum_{i=1}^{r}w_{i} \right) dw_{1} \cdots dw_{r} \\ \end{aligned}

となります。厳密には定積分と極限の交換は、被積分関数が一様収束することが十分条件となりますが、一様収束に関する証明は省略します。

次に \displaystyle\lim_{n \to \infty} \hat{m}_{X} を求めます。
\hat{m}_{X} については、 以下の式を満たします。

\begin{aligned} \dfrac{1}{\hat{m}_{X}} + \dfrac{1}{r}\sum_{i=1}^{r}\ln{x_{(i)}} -\dfrac{\sum_{i=1}^{r} x_{(i)}^{\hat{m}_{X}} \ln{x_{(i)}} +(n-r) x_{(r)}^{\hat{m}_{X}}\ln x_{(r)}} {\sum_{i=1}^{r} x_{(i)}^{\hat{m}_{X}} +(n-r) x_{(r)}^{\hat{m}_{X}} } =0 \end{aligned}

この式をまず、X_{(i)} \to Z_{(i)} に変数変換します。x_{(i)}=\eta z_{(i)}^{\frac{1}{m}} となることを利用して整理すると、

\begin{aligned} \dfrac{m}{\hat{m}_{Z}} + \dfrac{1}{r}\sum_{i=1}^{r}\ln{z_{(i)}} -\dfrac{\sum_{i=1}^{r} z_{(i)}^{\frac{\hat{m}_{Z}}{m}} \ln{z_{(i)}} +(n-r) z_{(r)}^{\frac{\hat{m}_{Z}}{m}}\ln z_{(r)}} {\sum_{i=1}^{r} z_{(i)}^{\frac{\hat{m}_{Z}}{m}} +(n-r) z_{(r)}^{\frac{\hat{m}_{Z}}{m}} } =0 \end{aligned}

となります。次に、Z_{(i)} \to W_{i} に変数変換します。z_{(i)}=\sum_{j=1}^{i} \dfrac{w_{j}}{n-j+1} となることから

\begin{aligned} \dfrac{m}{\hat{m}_{W}} + \dfrac{1}{r}\sum_{i=1}^{r}\ln \left( \sum_{j=1}^{i} \dfrac{w_{j}}{n-j+1} \right) -\dfrac{\sum_{i=1}^{r} \left( \sum_{j=1}^{i} \dfrac{w_{j}}{n-j+1} \right) ^{\frac{\hat{m}_{W}}{m}} \ln \left( \sum_{j=1}^{i} \dfrac{w_{j}}{n-j+1} \right) +(n-r) \left( \sum_{j=1}^{r} \dfrac{w_{j}}{n-j+1} \right)^{\frac{\hat{m}_{W}}{m}}\ln \left( \sum_{j=1}^{r} \dfrac{w_{j}}{n-j+1} \right)} {\sum_{i=1}^{r} \left( \sum_{j=1}^{i} \dfrac{w_{j}}{n-j+1} \right)^{\frac{\hat{m}_{W}}{m}} +(n-r) \left( \sum_{j=1}^{i} \dfrac{w_{j}}{n-j+1} \right)^{\frac{\hat{m}_{W}}{m}} } =0 \end{aligned}

と表せます。左辺の第3項については n が非常に大きいことから

\begin{aligned} \sum_{i=1}^{r} \left( \sum_{j=1}^{i} \dfrac{w_{j}}{n-j+1} \right) ^{\frac{\hat{m}_{W}}{m}} \ll (n-r) \left( \sum_{j=1}^{r} \dfrac{w_{j}}{n-j+1} \right)^{\frac{\hat{m}_{W}}{m}} \end{aligned}

となることを利用して、

\begin{aligned} &\dfrac{\sum_{i=1}^{r} \left( \sum_{j=1}^{i} \dfrac{w_{j}}{n-j+1} \right) ^{\frac{\hat{m}_{W}}{m}} \ln \left( \sum_{j=1}^{i} \dfrac{w_{j}}{n-j+1} \right) +(n-r) \left( \sum_{j=1}^{r} \dfrac{w_{j}}{n-j+1} \right)^{\frac{\hat{m}_{W}}{m}}\ln \left( \sum_{j=1}^{r} \dfrac{w_{j}}{n-j+1} \right)} {\sum_{i=1}^{r} \left( \sum_{j=1}^{i} \dfrac{w_{j}}{n-j+1} \right)^{\frac{\hat{m}_{W}}{m}} +(n-r) \left( \sum_{j=1}^{i} \dfrac{w_{j}}{n-j+1} \right)^{\frac{\hat{m}_{W}}{m}} } \\ &\approx \dfrac{ (n-r) \left( \sum_{j=1}^{r} \dfrac{w_{j}}{n-j+1} \right)^{\frac{\hat{m}_{W}}{m}}\ln \left( \sum_{j=1}^{r} \dfrac{w_{j}}{n-j+1} \right)} {(n-r) \left( \sum_{j=1}^{i} \dfrac{w_{j}}{n-j+1} \right)^{\frac{\hat{m}_{W}}{m}} } \\ &=\ln \left( \sum_{j=1}^{r} \dfrac{w_{j}}{n-j+1} \right) \end{aligned}

となります。よって式を変形すると、

\begin{aligned} \dfrac{m}{\hat{m}_{W}} &+ \dfrac{1}{r}\sum_{i=1}^{r}\ln \left( \sum_{j=1}^{i} \dfrac{w_{j}}{n-j+1} \right) -\ln \left( \sum_{j=1}^{r} \dfrac{w_{j}}{n-j+1} \right) \approx 0 \\ \dfrac{m}{\hat{m}_{W}} &\approx \ln \left( \sum_{j=1}^{r} \dfrac{w_{j}}{n-j+1} \right) -\dfrac{1}{r}\sum_{i=1}^{r} \ln \left( \sum_{j=1}^{i} \dfrac{w_{j}}{n-j+1} \right) \\ \hat{m}_{W} &\approx \dfrac{m}{\ln \left( \sum_{j=1}^{r} \dfrac{w_{j}}{n-j+1} \right) -\dfrac{1}{r}\sum_{i=1}^{r} \ln \left( \sum_{j=1}^{i} \dfrac{w_{j}}{n-j+1} \right) } \\ &= \dfrac{mr}{r \ln \left( \sum_{j=1}^{r} \dfrac{w_{j}}{n-j+1} \right) -\sum_{i=1}^{r} \ln \left( \sum_{j=1}^{i} \dfrac{w_{j}}{n-j+1} \right) } \\ &= \dfrac{mr}{ \ln \left\{ \dfrac{\left( \sum_{j=1}^{r} \dfrac{w_{j}}{n-j+1} \right)^r} { \prod_{i=1}^{r} \left( \sum_{j=1}^{i} \dfrac{w_{j}}{n-j+1} \right) } \right\}} \\ &= \dfrac{mr}{ \ln \left\{ \left( \dfrac{\sum_{j=1}^{r} \dfrac{w_{j}}{n-j+1}}{\dfrac{w_{1}}{n}} \right) \cdots \left( \dfrac{\sum_{j=1}^{r} \dfrac{w_{j}}{n-j+1}}{\sum_{j=1}^{r-1} \dfrac{w_{j}}{n-j+1}} \right) \right\}} \\ \end{aligned}

となります。
よって、\hat{m}_{W} について n\to \infty の極限をとると、以下の通り \displaystyle \lim_{n \to \infty} \hat{m}_{W} を求めることができます。

\begin{aligned} \lim_{n \to \infty} \hat{m}_{W} &= \lim_{n \to \infty}\dfrac{mr}{ \ln \left\{ \left( \dfrac{\sum_{j=1}^{r} \dfrac{nw_{j}}{n-j+1}}{w_{1}} \right) \cdots \left( \dfrac{\sum_{j=1}^{r} \dfrac{nw_{j}}{n-j+1}}{\sum_{j=1}^{r-1} \dfrac{nw_{j}}{n-j+1}} \right) \right\}} \\ &= \dfrac{mr}{ \ln \left( \dfrac{\sum_{j=1}^{r} w_{j} }{ w_{1}} \cdots \dfrac{\sum_{j=1}^{r} w_{j} }{\sum_{j=1}^{r-1} w_{j }} \right) } \\ \end{aligned}

以上より、目標の \displaystyle \lim_{n\to \infty} E[\hat{m}_{W}] については、

\begin{aligned} \lim_{n\to \infty }E[\hat{m_{X}}] &=\int_{0}^{\infty} \cdots \int_{0}^{\infty} \dfrac{mr}{ \ln \left( \dfrac{\sum_{i=1}^{r} w_{i} }{ w_{i}} \cdots \dfrac{\sum_{i=1}^{r} w_{i} }{\sum_{i=1}^{r-1} w_{i}} \right) } \exp\left( -\sum_{i=1}^{r}w_{i} \right) dw_{1} \cdots dw_{r} \\ \end{aligned}

と表せます。W_{i} \ (i=1,\ ...\ ,r) は互いに独立に母数 1 の指数分布に従います。

2.4. W1, ... , Wr から U(1), ... , U(r-1),V への変数変換

次に、以下のように W_{1},\ ...\ ,W_{r} \to U_{(1)},\ ... \ ,U_{(r-1)}, V の変数変換します。

\begin{aligned} &V=W_{1}+ \cdots +W_{r} \\ &U_{(i)} = \dfrac{\sum_{j=1}^{i} W_{j}}{V}\hspace{10pt} (i=1,\ ...\ ,r-1) \\ \end{aligned}

W_{i}U_{(i)},V の形に変形すると、U_{(i)}V =W_{1}+ \cdots + W_{i} =U_{(i-1)}V+W_{i} から、

\left\{\begin{aligned} &W_{1} = VU_{(1)} \\ &W_{i} = V( U_{(i)} -U_{(i-1)} ) \hspace{10pt} (i=2,\ ...\ ,r-1) \\ &W_{r} = V(1-U_{(r-1)}) \end{aligned}\right.

となります。よってヤコビアン J

\begin{aligned} J&=\begin{vmatrix} \dfrac{\partial w_{1}}{\partial u_{(1)}} & \dfrac{\partial w_{2}}{\partial u_{(1)}} & \cdots & \dfrac{\partial w_{r-1}}{\partial u_{(1)}} &\dfrac{\partial w_{r}}{\partial u_{(1)}} \\\\ \dfrac{\partial w_{1}}{\partial u_{(2)}} & \dfrac{\partial w_{2}}{\partial u_{(2)}} & \cdots & \dfrac{\partial w_{r-1}}{\partial u_{(2)}} &\dfrac{\partial w_{r}}{\partial u_{(2)}} \\\\ \vdots & \vdots & \ddots & \vdots & \vdots \\\\ \dfrac{\partial w_{1}}{\partial u_{(r-1)}} & \dfrac{\partial w_{2}}{\partial u_{(r-1)}} & \cdots & \dfrac{\partial w_{r-1}}{\partial u_{(r-1)}} &\dfrac{\partial w_{r}}{\partial u_{(r-1)}} \\\\ \dfrac{\partial w_{1}}{\partial v} & \dfrac{\partial w_{2}}{\partial v} & \cdots & \dfrac{\partial w_{r-1}}{\partial v} &\dfrac{\partial w_{r}}{\partial v} \\\\ \end{vmatrix} \\ J&=\begin{vmatrix} v & -v & \cdots & 0 & 0 \\ 0 & v & \cdots & 0 & 0 \\ \vdots & \vdots & \ddots & \vdots & \vdots \\ 0 & 0 & \cdots & v & -v \\ u_{(1)} & u_{(2)} - u_{(1)} & \cdots & u_{(r-1)} - u_{(r-2)} & 1-u_{(r-1)} \\\\ \end{vmatrix} \\ &=\begin{vmatrix} v & 0 & \cdots & 0 & 0 \\ 0 & v & \cdots & 0 & 0 \\ \vdots & \vdots & \ddots & \vdots & \vdots \\ 0 & 0 & \cdots & v & -v \\ u_{(1)} & u_{(2)} & \cdots & u_{(r-1)} - u_{(r-2)} & 1-u_{(r-1)} \\ \end{vmatrix} \\ &= \cdots \\ &=\begin{vmatrix} v & 0 & \cdots & 0 & 0 \\ 0 & v & \cdots & 0 & 0 \\ \vdots & \vdots & \ddots & \vdots & \vdots \\ 0 & 0 & \cdots & v & 0 \\ u_{(1)} & u_{(2)} & \cdots & u_{(r-1)} & 1 \\ \end{vmatrix} \\ &= v^{r-1} \\ \end{aligned}

なお、u_{(i)} の積分区間については w_{i} \gt 0 より

\begin{aligned} 0 \lt u_{(1)} \lt \cdots \lt u_{(r-1)} \lt 1 \end{aligned}

と表せます。よって求める期待値 \displaystyle \lim_{n \to \infty} E[\hat{m}_{X}]

\begin{aligned} \lim_{n \rightarrow \infty}E[\hat{m}_{X}] &=\int_{0}^{1} \cdots \int_{0}^{u_{(2)}} \int_{0}^{\infty} \dfrac{mr}{\ln \left( \dfrac{1}{u_{(1)} \cdots u_{(r-1)}} \right) } \exp\left(-v \right) v^{r-1} dv du_{(1)} \cdots du_{(r-1)} \\ &= \int_{0}^{1} \cdots \int_{0}^{u_{(2)}} -\dfrac{ mr\cdot (r-1)!}{ \ln \left( u_{(1)} \cdots u_{(r-1)} \right) } \int_{0}^{\infty} \dfrac{v^{r-1}}{(r-1)!} \exp\left(-v \right) dv du_{(1)} \cdots du_{(r-1)} \\ &= mr \int_{0}^{1} \cdots \int_{0}^{u_{(2)}} -\dfrac{1}{ \ln \left( u_{(1)} \cdots u_{(r-1)} \right) } (r-1)!\ du_{(1)} \cdots du_{(r-1)} \\ \end{aligned}

と変形できます。

\displaystyle \int_{0}^{\infty} \dfrac{v^{r-1}}{(r-1)!} \exp\left(-v \right) dv については、Gamma(r,1) に従う分布において、すべての定義域で積分することから 1 となります。

2.5. U(1), ... , U(r-1) から T2, ... , Tr-1,S への変数変換

次は以下のように、U_{(1)},\ ...\ ,U_{(r-1)} \to T_{2},\ ...\ ,T_{r-1},S の変数変換をします。

\begin{aligned} &T_{i} = -\ln U_{(i)}\ (i=2,\ ...\ ,r-1) \\ S&= - \ln(U_{(1)} \cdots U_{(r-1)}) \\ \end{aligned}

U_{(1)},\ ...\ ,U_{(r-1)}T_{2},\ ...\ ,T_{r-1},S で表すと、

\begin{aligned} U_{(i)} &= e^{-T_{i}} (i=2,\ ...\ ,r-1) \\ &=-\ln U_{(1)} -\ln U_{(2)} \cdots - \ln U_{(r-1)} \\ &= -\ln U_{(1)} + T_{2}+\cdots +T_{r-1} \\ &\ln U_{(1)} = T_{2}+ \cdots +T_{r-1} -S\\ & U_{(1)} =e^{T_{2}+\cdots+T_{r-1}-S} \\ \end{aligned}

となります。ヤコビアン J については、

\begin{aligned} J&=\begin{vmatrix} \dfrac{\partial u_{(1)}}{\partial s} & \dfrac{\partial u_{(2)}}{\partial s} & \cdots & \dfrac{\partial u_{(r-2)}}{\partial s} &\dfrac{\partial u_{(r-1)}}{\partial s} \\\\ \dfrac{\partial u_{(1)}}{\partial t_{2}} & \dfrac{\partial u_{(2)}}{\partial t_{2}} & \cdots & \dfrac{\partial u_{(r-2)}}{\partial t_{2}} &\dfrac{\partial u_{(r-1)}}{\partial t_{2}} \\\\ \vdots & \vdots & \ddots & \vdots & \vdots \\\\ \dfrac{\partial u_{(1)}}{\partial t_{r-2}} & \dfrac{\partial u_{(2)}}{\partial _{r-2}} & \cdots & \dfrac{\partial u_{(r-2)}}{\partial t_{r-2}} &\dfrac{\partial u_{(r-1)}}{\partial t_{r-2}} \\\\ \dfrac{\partial u_{(1)}}{\partial t_{r-1}} & \dfrac{\partial u_{(2)}}{\partial t_{r-1}} & \cdots & \dfrac{\partial u_{(r-2)}}{\partial t_{r-1}} &\dfrac{\partial u_{(r-1)}}{\partial t_{r-1}} \\\\ \end{vmatrix}\\ &=\begin{vmatrix} -e^{t_{2}+\cdots+t_{r-1}-s} & 0 & \cdots & 0 & 0 \\ e^{t_{2}+\cdots+t_{r-2}-s} & -e^{-t_{2}} & \cdots & 0 & 0 \\ \vdots & \vdots & \ddots & \vdots & \vdots \\ e^{t_{2}+\cdots+t_{r-2}-s} & 0 & \cdots & -e^{-t_{r-2}} & 0 \\ e^{t_{2}+\cdots+t_{r-2}-s} & 0 & \cdots & 0 & -e^{-t_{r-1}} \\ \end{vmatrix} \\ &=\left| (-1)^{r-1} e^{-s} \right| \\ &=e^{-s} \end{aligned}

積分区間については 0 \lt u_{(1)} \lt \cdots \lt u_{(r-1)} \lt 1 より

\begin{aligned} &0 \lt e^{t_{2}+\cdots +t_{r-1}-s} \lt e^{-t_{2}} \lt \cdots \lt e^{-t_{r-1}} \lt 1 \\ & 0 \lt t_{r-1} \lt \cdots \lt t_{2} \lt s-t_{2}-\cdots -t_{r-1} \\ & 0 \lt t_{r-1}\lt \cdots \lt t_{4} \lt t_{3} \lt t_{2} \lt \dfrac{s-t_{3}-\cdots - t_{r-1}}{2} \\ \end{aligned}

となります。t_{3} の積分区間については以下のようになります。

\begin{aligned} t_{4} \lt t_{3} \lt \dfrac{s-t_{3}-\cdots - t_{r-1}}{2} \\ 2t_{3} \lt s-t_{3}-\cdots - t_{r-1} \\ 3t_{3} \lt s-t_{4}-\cdots - t_{r-1} \\ t_{4} \lt t_{3} \lt \dfrac{s-t_{4}-\cdots - t_{r-1}}{3} \\ \end{aligned}

よって、t_{i} (i=3,\ ...\ ,r-2) の積分区間については

\begin{aligned} t_{i+1} \lt t_{i} \lt \dfrac{s-t_{i}-\cdots - t_{r-1}}{i-1} \\ (i-1)t_{i} \lt s-t_{i}-\cdots - t_{r-1} \\ i t_{i} \lt s-t_{i+1}-\cdots - t_{r-1} \\ t_{i+1} \lt t_{i} \lt \dfrac{s-t_{i+1}-\cdots - t_{r-1}}{i} \\ \end{aligned}

と表せます。以上より、求める期待値 \displaystyle \lim_{n \to \infty} E[\hat{m}_{X}]

\begin{aligned} \lim_{n \to \infty} E[\hat{m}_{X}] &=mr\int_{0}^{\infty} \dfrac{e^{-s}}{s} (r-1)! \int_{0}^{\frac{s}{r-1}} \cdots \int_{t_{3}}^{\frac{s-t_{3}-\cdots - t_{r-1}}{2}} dt_{2} \cdots dt_{r-1} ds \\ &=mr\int_{0}^{\infty} \dfrac{e^{-s}}{s} (r-1)! \int_{0}^{\frac{s}{r-1}} \cdots \int_{t_{4}}^{\frac{s-t_{4}-\cdots - t_{r-1}}{3}} \dfrac{s-3t_{3} -t_{4} - \cdots -t_{r-1}}{2} dt_{3} \cdots dt_{r-1} ds \\ &=mr\int_{0}^{\infty} \dfrac{e^{-s}}{s} (r-1)! \int_{0}^{\frac{s}{r-1}} \cdots \int_{t_{5}}^{\frac{s-t_{5}-\cdots - t_{r-1}}{4}} -\left[ \dfrac{1}{3 \cdot 2 \cdot 2! 1!} (s-3t_{3} -t_{4} - \cdots -t_{r-1} )^{2} \right]_{t_{4}}^{\frac{s-t_{4}-\cdots - t_{r-1}}{3}} dt_{4} \cdots dt_{r-1} ds \\ &=mr\int_{0}^{\infty} \dfrac{e^{-s}}{s} (r-1)! \int_{0}^{\frac{s}{r-1}} \cdots \int_{t_{5}}^{\frac{s-t_{5}-\cdots - t_{r-1}}{4}} \dfrac{1}{3! 2!} (s-4t_{4} -t_{5} - \cdots -t_{r-1} )^{2} dt_{4} \cdots dt_{r-1} ds \\ &=\cdots \\ &=mr\int_{0}^{\infty} \dfrac{e^{-s}}{s} (r-1)! \int_{0}^{\frac{s}{r-1}} \dfrac{1}{(r-2)! (r-3)!} (s-(r-1)t_{r-1})^{r-3} dt_{r-1} ds \\ &=mr\int_{0}^{\infty} \dfrac{e^{-s}}{s} (r-1)! \left[ \dfrac{1}{ (r-1)(r-2) (r-2)!(r-3)!} (s-(r-1)t_{r-1})^{r-2} \right]_{0}^{\frac{s}{r-1}} dt_{r-1} ds \\ &=mr\int_{0}^{\infty} \dfrac{e^{-s}}{s} (r-1)! \dfrac{s^{r-2}}{ (r-1)! (r-2)! } ds \\ &=mr\dfrac{1}{(r-2)!} \int_{0}^{\infty} e^{-s}s^{r-3} ds \\ &=mr\dfrac{(r-3)!}{(r-2)!} = \dfrac{mr}{r-2} \end{aligned}

となり、以下の式が導出されます。

\begin{aligned} \lim_{n \to \infty} E[\hat{m}_{X}] = \dfrac{r}{r-2}m \end{aligned}

k_{r,n}=\dfrac{\hat{m}_{X}}{m} とすると、

\begin{aligned} \lim_{n\to \infty}E[k_{r,n}]=\dfrac{r}{r-2} \end{aligned}

が成り立ちます。

Discussion