ライントランセクト調査による個体数の区間推定4
今回は$ \bar{s} $の分散について説明します。
まず、何も調整を行わずに、単純に得られたデータの群れサイズを平均する場合です。この場合は、
\bar{s} = \frac{1}{n}\sum_{i=1}^{n} s_i
とするだけです。siはi番目の群れサイズを表します。分散は、
\begin{aligned} Var[\bar{s}] = \frac{1}{n^2}\sum_{i=1}^{n} Var[s_i] & = \frac{1}{n^2}・n・\frac{1}{n-1} \sum_{i=1}^{n} \left( s_i - \bar{s} \right)^2 \\ & = \frac{1}{n(n-1)} \sum_{i=1}^n \left(s_i - \bar{s} \right)^2 \end{aligned}
となります。
しかしこれには問題があります。通常ライントランセクト調査では、トラックラインから遠い場所だと、群れサイズが大きい群ればかり発見され、群れサイズが小さい群れのほうが見逃しやすいという事態が発生するのです。
このような場合、単純にデータの群れサイズを平均しただけでは群れサイズを過小推定する恐れがあるため、群れサイズの調整を行うことが多いです。
ではどのように調整すればよいのでしょうか。実際には、線形回帰分析により平均群れサイズを推定する手法がしばしば用いられます(Buckland et al., 1993)。モデル式は以下の通りです。
z = \log(s) = a+bg(x)+\epsilon, \hspace{8mm} \epsilon ~ \mathcal{N}(0,\sigma^2)
つまり、データ上では$ z=\log(s)とg(x) $は一次関数の関係になっていて、誤差項は正規分布に従う、ということを仮定するのです。
ここからどのように群れサイズ平均を推定するかというと、すべての個体が発見される場合、つまりg(x)=1の地点におけるsの確率分布が、真の群れサイズの確率分布と同じになると考えるのです。したがって、$ g(x)=1 $におけるsの確率分布を考えれば良いということになります。
大半の人は回帰分析を大学1年で習うので基礎的な部分は省略しますが、パラメータaとbの推定値は以下のようになります。
\hat{a} = \bar{z} - b \bar{g}(x) \tag{1}
\hat{b} = \frac{S_{xy}}{S_{xx}}
ただし、
S_{xy} = \frac{1}{n} \sum_{i=1}^{n} (g(x_i) - \bar{g}(x))(z_i - \bar{z})
S_{xx} = \frac{1}{n} \sum_{i=1}^{n} (g(x_i) - \bar{g}(x))^2
です。
モデル式より、g(x)=1の時のときのzの$ z_0 $とすると、
z_0 = a+b+\epsilon
(1)のaの推定値を代入すると、
z_0 = \bar{z} - b \bar{g}(x) +b+\epsilon
\hspace{8mm} = \bar{z}+b(1-\bar{g}(x)) +\epsilon
よって、
\begin{aligned}Var[z_0] = & Var[\bar{z}] + (1-\bar{g}(x))^{2}Var[b] + Var[\epsilon] + \\ & 2(1-\bar{g}(x))Cov[\bar{z},b] + 2(1-\bar{g}(x))Cov[b,\epsilon] + 2Cov[\bar{z},\epsilon] \end{aligned} \tag{2}
となります。ただしここでは、$ \epsilon $はデータと独立なため、$ C ov[b, \epsilon], C ov[\bar{z}, \epsilon] $はともに0になります。
回帰分析に関する本はたくさんあるので詳しい説明はそちらに譲りますが、各分散、共分散の値は以下のようになります。
Var[\bar{z}] = \frac{\sigma^2}{n}
Var[b] = \frac{\sigma^2}{S_{xx}}
Var[\epsilon] = \sigma^2
Cov[\bar{z},b] = 0
となります。ただし、$ σ^2 $は残差分散で、その不偏推定量は
\hat{\sigma^2} = \frac{1}{n-2} \sum_{i=1}^{n} (z_i - a+bg(x_i))^2
と求められます。
これらの値を式(2)に代入して整理すると、
Var[z_0] = \sigma^2 \left(1 + \frac{1}{n} + \frac{(1-\bar{g}(x))^{2}}{S_{xx}} \right)
となります。
以上より、zは期待値$ a+b $、分散$ Var[z_0] $の正規分布に従うことが分かりました。しかし今求めたいのは、$ s_0=\exp(z_0) $の確率分布です。
(正規分布に従う確率変数を指数変換した時、その分布は対数正規分布と呼ばれるものになります。)
期待値と分散を考えてみます。これらは、正規分布の確率密度関数まで戻り、期待値と分散の定義に従って積分で求めることが出来ます。
E[s_0] = E[\exp(z_0)] = \int_{-∞}^{∞}\frac{\exp(z)}{\sqrt{2 \pi Var[z_0]}}\exp \left(-\frac{(z-(a+b))^2}{2Var[z_0]} \right) dx
途中計算は面倒なので省略します。exp()の中身を平方完成して整理していくと、以下の式が導けます。
E[s_0] = \exp \left(a+b+\frac{Var[z_0]}{2} \right)
同様に$ E[s_0^2] $は
E[{s_0}^2] = \exp(2(a+b)+2Var[z_0])
となります。故に、
Var[s_0] = E[{s_0}^2] = E[s_0]^2 = \exp(2(a+b)+Var[z_0])・(\exp(Var[z_0])-1)
となります。これが、$ g(x)=1 $となる点における群れサイズsの確率分布の期待値と分散です。
これをもとに群れサイズの平均について考えれば良いわけです。平均の期待値は$ E[s_0] $に一致します。分散は、上で求めた$ s_0 $の分散をnで割れば良いだけです。