非对称损失函数下对数伽玛分布尺度参数的Bayes估计

Bayes Estimation of Scale Parameter of Log Gamma Distribution under Asymmetric Loss Function

摘要: 损失函数和先验分布的选取在Bayes估计问题中起重要作用。已有众多学者在不同损失函数，不同先验分布下研究了对数伽玛分布尺度参数的Bayes估计问题。本文基于Quasi先验分布分别在Linex损失与熵损失函数下，研究了对数伽玛分布尺度参数的Bayes估计。并给出了其在熵损失函数下的E-Bayes估计和多层Bayes估计。最后通过数值模拟方法对各种估计结果的优良性进行比较分析，结果表明，先验分布为Quasi分布时，尺度参数θ的Linex损失函数下的Bayes估计较熵损失函数下的稳健性更好，更接近真值。

Abstract: The selection of loss function and prior distribution plays an important role in Bayes estimation. Many scholars have studied the Bayes estimation of log-gamma distribution scale parameter under different loss function and different prior distribution. Based on Quasi-prior distribution and Linex loss function and entropy loss function respectively, the Bayes of log-gamma distribution scale parameter is studied. And its E-Bayes estimation and multilayer Bayes estimation under entropy loss function are given. Finally, the numerical simulation method is used to compare the robustness of various estimation results. The results show that when the prior distribution is Quasi, the Bayes estimation under the Linex loss function of the scale parameter theta is more stable and accurate than that under the entropy loss function.

1. 引言

在统计决策理论中，参数估计的优劣很大程度上依赖于损失函数形式的选择。Linex损失函数自Varian [1] 提出以来，已被众多学者应用于参数估计中。王琪等 [2] 研究了Linex和熵损失这两类非对称损失函数下逆指数分布参数的Bayes估计。王成元等 [3] 对对数伽玛分布尺度参数在Linex损失与复合Linex损失函数下的Bayes估计进行了比较。金秀岩 [4] 研究了复合MLinex对称损失函数下对数伽玛分布参数的Bayes估计。对数伽玛分布常应用于风险理论 [5]、油气资源评价 [6] 等领域。

设随机变量服从对数伽玛分布，其概率密度函数为

$f\left(x;\theta ,\alpha \right)=\frac{{\theta }^{\alpha }{\left(\ln x\right)}^{\alpha -1}}{\Gamma \left(\alpha \right)}{x}^{-\left(\theta +1\right)},x>0,\alpha >0,\theta >0$ (1)

设 $X_1,X_2,\cdots ,X_n$ 为来自对数伽玛模型的一个简单样本，则对其n个独立观测值 $x_1,\cdots ,x_n$ 得到似然函数为

$\begin{array}{c}f\left(x_1,x_2,\cdots ,x_n;\theta \right)=\frac{{\theta }^{n\alpha }}{{\left[\Gamma \left(\alpha \right)\right]}^n}{\displaystyle \underset{i=1}{\overset{n}{\prod }}{\left(\ln x_i\right)}^{\alpha -1}}{\left(x_1{x}_2\cdots x_n\right)}^{-\left(\theta +1\right)}\\ =\frac{{\theta }^{n\alpha }}{{\left[\Gamma \left(\alpha \right)\right]}^n}{\displaystyle \underset{i=1}{\overset{n}{\prod }}{\left(\ln x_i\right)}^{\alpha -1}}{\text{e}}^{-\left(\theta +1\right){\displaystyle \underset{i=1}{\overset{n}{\sum }}\ln x_i}}\propto {\theta }^{n\alpha }{\text{e}}^{-\theta t}\end{array}$ (2)

其中 $t={\displaystyle \underset{i=1}{\overset{n}{\sum }}\ln x_i}$。

设参数θ的先验分布为无信息Quasi先验分布，相应的概率密度函数为

$\pi \left(\theta \right)\propto \frac{1}{{\theta }^d},\theta >0,d>0$ (3)

当 $d=0$ 时， $\pi \left(\theta \right)\propto 1$ 为离散先验分布；当 $d=1$ 时， $\pi \left(\theta \right)\propto \frac{1}{\theta }$ 为无信息先验分布。

定义1 [1] 定义Linex损失函数如下式

$L_a\left(\theta ,\delta \right)={\text{e}}^{a\left(\delta -\theta \right)}-a\left(\delta -\theta \right)-1,\left(a>0\right)$ (4)

其中 $\delta$ 为参数θ的估计，a是该损失函数的尺度参数， $a\ne 0$。

事实上对 $\delta$ 求偏导得 $L^{\prime }\left(\theta ,\delta \right)=a{\text{e}}^{a\left(\delta -\theta \right)}-a$，任取 $0<{\delta }_1<{\delta }_2$ 则有 $L^{\prime }\left(\theta ,{\delta }_1\right)-L^{\prime }\left(\theta ,{\delta }_2\right)=a\left({\text{e}}^{a\left({\delta }_1-\theta \right)}-{\text{e}}^{a\left({\delta }_2-\theta \right)}\right)<0$，所以Linex损失函数关于 $\delta$ 是严格凸函数。

定义2 [7] 定义熵损失函数如下式

$L\left(\theta ,\delta \right)=\frac{\delta }{\theta }-\ln \frac{\delta }{\theta }-1$ (5)

事实上对 $\delta$ 求偏导得 $L^{\prime }\left(\theta ,\delta \right)=\frac{1}{\theta }-\frac{1}{\delta }$，任取 $0<{\delta }_1<{\delta }_2$ 则有 $L^{\prime }\left(\theta ,{\delta }_1\right)-L^{\prime }\left(\theta ,{\delta }_2\right)=\frac{{\delta }_1-{\delta }_2}{{\delta }_1{\delta }_2}<0$，所以熵损失函数对于 $\delta$ 是严格凸函数。

定义3 [8] 对 $\left(\beta ,\lambda \right)\in D$，若 ${\delta }_B\left(\beta ,\lambda \right)$ 是连续的，则称

${\delta }_{EB}={\displaystyle {\iint }_D{\delta }_B}\left(\beta ,\lambda \right)\pi \left(\beta ,\lambda \right)\text{d}\beta \text{d}\lambda$ (6)

是参数θ的E-Bayes估计。其中 ${\displaystyle {\iint }_D{\delta }_B}\left(\beta ,\lambda \right)\pi \left(\beta ,\lambda \right)\text{d}\beta \text{d}\lambda$ 是存在的，D为超参数 $\beta$ 和 $\lambda$ 取值的集合， $\pi \left(\beta ,\lambda \right)$ 是 $\beta$ 和 $\lambda$ 在集合D上的密度函数， ${\delta }_B\left(\beta ,\lambda \right)$ 为 $\theta$ 的Bayes估计。

引理1 [9] 在Linex损失函数(4)下，对任何先验分布 $\pi \left(\theta \right)$， $\theta$ 唯一的Bayes估计为

${\delta }_B\left(x\right)=-\frac{1}{a}\ln \left(E\left({\text{e}}^{-a\theta }|X\right)\right)$ (7)

引理2 [2] 在熵损失函数(5)下，对任何先验分布 $\pi \left(\theta \right)$，θ唯一的Bayes估计为

${\delta }_B\left(x\right)={\left[E\left({\theta }^{-1}|X\right)\right]}^{-1}$ (8)

关于对数伽玛分布尺度参数的先验分布为Quasi分布时，在非对称损失函数下的Bayes估计文献研究较少，为此，在Linex损失与熵损失这两种非对称损失函数下研究了对数伽玛分布尺度参数的Bayes估计，最后得出：当先验分布同为Quasi先验分布时，Linex损失函数下的Bayes估计的稳健性较熵损失函数下的Bayes估计稳健性更好，更接近真值。本文共有三部分，第一部分主要介绍研究背景；第二部分为主要结果，在Linex损失函数，熵损失函数下研究了参数θ的Bayes估计以及当先验分布改为伽玛分布时，推导出了尺度参数θ在熵损失函数下的E-Bayes估计和多层Bayes估计；最后一个部分通过R软件进行数值模拟，比较了各估计值的精确性和稳健性。

2. 主要结果

2.1. 参数θ的Bayes估计

在损失函数(4)，(5)下，结合引理1、引理2讨论对数伽玛分布的尺度参数θ的Bayes估计问题。

引理3 设随机变量X服从对数伽玛分布(1)，参数θ的先验分布 $\pi \left(\theta \right)$ 为Quasi先验分布，

$X_1,X_2,\cdots ,X_n$ 是服从对数伽玛分布的一个简单样本，记 $X=\left(X_1,\cdot \cdot \cdot ,X_n\right)$， $T={\displaystyle \underset{i=1}{\overset{n}{\sum }}\ln x_i}$ 则θ的后验分布服从

$Gamma\left(n\alpha -d+1,T\right)$。

证明 由(2)式可知 $X_1,X_2,\cdots ,X_n$ 的n个独立观察值 $x_1,\cdots ,x_n$ 的似然函数。再由(3)式及Bayes定理，参数θ的后验密度函数为

$f\left(\theta |x\right)\propto f\left(x_1,\cdots ,x_n;\theta \right)\pi \left(\theta \right)\propto {\theta }^{n\alpha -d}{\text{e}}^{-\theta t}$ (9)

于是 $\theta |x~Gamma\left(n\alpha -d+1,t\right)$，所以 $\theta |X~Gamma\left(n\alpha -d+1,T\right)$。

引理4 [10] 在给定的Bayes统计决策问题中，假如对给定的先验分布 $\pi \left(\theta \right)$，θ的Bayes估计 ${\delta }_B$ 是唯一的，则它也是容许的。

定理1在Linex非对称损失函数下，假设对数伽玛分(1)的尺度参数θ的先验分布为Quasi先验分布，则尺度参数θ的Bayes估计为

${\delta }_B\left(x\right)=\frac{n\alpha -d+1}{a}\ln \frac{a}{T}$ (10)

证明 由引理3知参数θ的后验密度函数为

$H\left(\theta |x\right)=\frac{T^{n\alpha -d+1}}{\Gamma \left(n\alpha -d+1\right)}{\theta }^{n\alpha -d}{\text{e}}^{-T\theta }$ (11)

又

$\begin{array}{c}E\left({\text{e}}^{-a\theta }|X\right)={\displaystyle {\int }_0^{\infty }{\text{e}}^{-a\theta }}\frac{T^{n\alpha -d+1}}{\Gamma \left(n\alpha -d+1\right)}{\theta }^{n\alpha -d}{\text{e}}^{-T\theta }\text{d}\theta \\ ={\displaystyle {\int }_0^{\infty }\frac{T^{n\alpha -d+1}}{\Gamma \left(n\alpha -d+1\right)}{\theta }^{n\alpha -d}{\text{e}}^{-\left(a+T\right)\theta }\text{d}\theta }\\ =\frac{T^{n\alpha -d+1}}{{\left(a+T\right)}^{n\alpha -d+1}}{\displaystyle {\int }_0^{\infty }\frac{{\left(a+T\right)}^{n\alpha -d+1}}{\Gamma \left(n\alpha -d+1\right)}{\theta }^{n\alpha -d}{\text{e}}^{-\left(a+T\right)\theta }\text{d}\theta }\\ ={\left(\frac{T}{a+T}\right)}^{n\alpha -d+1}\end{array}$ (12)

则由引理1得尺度参数θ的Bayes估计为

${\delta }_B\left(x\right)=E\left({\text{e}}^{-a\theta }|X\right)=-\frac{1}{a}\ln {\left(\frac{T}{a+T}\right)}^{n\alpha -d+1}=\frac{n\alpha -d+1}{a}\ln \frac{T+a}{T}$ (13)

定理2在熵损失函数下，假设对数伽玛分布(1)的尺度参数θ的先验分布为Quasi先验分布，则尺度参数θ的Bayes估计为

${\delta }_B\left(x\right)=\frac{n\alpha -d}{T}$ (14)

运用引理2证明，证明过程类似定理1。

定理3在熵损失函数下，对数伽玛分布(1)的尺度参数θ的先验分布为伽玛分布 $Gamma\left(\beta ,\lambda \right)$，则尺度参数θ的Bayes估计为

${\delta }_B\left(x\right)=\frac{n\alpha +\beta -1}{T+\lambda }$ (15)

2.2. 参数θ的E-Bayes估计

在定理3的前提下讨论尺度参数θ的E-Bayes估计。先验分布为伽玛分布时，根据文献 [11] 的减函数法， $\beta$ 和 $\lambda$ 的选取应使先验分布 $\pi \left(\theta ;\beta ,\lambda \right)=\frac{{\lambda }^{\beta }}{\Gamma \left(\beta \right)}{\theta }^{\beta -1}\exp \left\{-\lambda \theta \right\}$ 为 $\theta$ 的单调减函数。

因为 $\pi \left(\theta ;\beta ,\lambda \right)$ 对θ求导数为 $\frac{\text{d}\pi \left(\theta ;\beta ,\lambda \right)}{\text{d}\theta }=\frac{{\lambda }^{\beta }{\theta }^{\beta -2}{\text{e}}^{-\lambda \theta }}{\Gamma \left(\beta \right)}\left[\left(\beta -1\right)-\lambda \theta \right]$。由伽玛分布的定义知， $\beta >0$， $\lambda >0$， $\theta >0$。所以，当 $0<\beta <1$， $\lambda >0$ 时， $\frac{\text{d}\pi \left(\theta ;\beta ,\lambda \right)}{\text{d}\theta }<0$，即 $\pi \left(\theta ;\beta ,\lambda \right)$ 为θ的单调减函数。然后由定义3，超参数 $\beta ,\lambda$ 的先验分布可以采取均匀分布的形式，即

$\pi \left(\beta \right)=U\left(0,1\right)$， $\pi \left(\lambda \right)=U\left(0,c\right)$，其中 $c>0$ 为 $\lambda$ 的上界。

定理4在熵损失函数下，对于对数伽玛分布(1)，若尺度参数θ的先验分布 $Gamma\left(\beta ,\lambda \right)=\frac{{\lambda }^{\beta }}{\Gamma \left(\beta \right)}{\theta }^{\beta -1}{\text{e}}^{-\lambda \theta }$ 的超参数 $\beta ,\lambda$ 的先验分布为D上的均匀分布，则尺度参数θ的E-Bayes估计为

${\delta }_{EB}\left(x\right)=\frac{1}{c}\left(n\alpha -\frac{1}{2}\right)\ln \left(\frac{T+c}{T}\right)$ (16)

证明 由假设知超参数 $\beta ,\lambda$ 的先验分布为 $\pi \left(\beta ,\lambda \right)=\frac{1}{c}$，其中 ${\delta }_B\left(x\right)$ 是由定理3得知，即是在熵损失函数下，对数伽玛分布(1)的尺度参数θ的先验分布为伽玛分布时，尺度参数θ的Bayes估计，其为 ${\delta }_B\left(x\right)=\frac{n\alpha +\beta -1}{T+\lambda }$。根据定义3得

$\begin{array}{c}{\delta }_{EB}\left(x\right)={\displaystyle {\iint }_D{\delta }_B\left(x\right)\pi \left(\beta ,\lambda \right)\text{d}\beta \text{d}\lambda }\\ ={\displaystyle {\int }_0^1{\displaystyle {\int }_0^c\frac{1}{c}\frac{n\alpha +\beta -1}{\lambda +T}\text{d}\lambda \text{d}\beta }}\\ =\frac{1}{c}{\displaystyle {\int }_0^1\left(n\alpha +\beta -1\right){\displaystyle {\int }_0^c\frac{1}{\lambda +T}\text{d}\lambda \text{d}\beta }}\\ =\frac{1}{c}\ln \left(\frac{T+c}{T}\right){\displaystyle {\int }_0^1\left(n\alpha +\beta -1\right)\text{d}\beta }\\ =\frac{1}{c}\left(n\alpha -\frac{1}{2}\right)\ln \left(\frac{T+c}{T}\right)\end{array}$

2.3. 参数θ的多层Bayes估计

最后根据上面的讨论推导出尺度参数θ的多层Bayes估计。

若参数θ的先验分布为 $Gamma\left(\beta ,\lambda \right)$，将 $\left(0,1\right)$， $\left(0,c\right)$ 上的均匀分布作为超参数 $\beta$ 和 $\lambda$ 的先验分布，则θ的多层先验密度函数为

$\pi \left(\theta \right)={\displaystyle {\int }_0^1{\displaystyle {\int }_0^c\pi \left(\theta ;\beta ,\lambda \right)\pi \left(\beta \right)\pi \left(\lambda \right)\text{d}\lambda \text{d}\beta }}=\frac{1}{c}{\displaystyle {\int }_0^1{\displaystyle {\int }_0^c\frac{{\lambda }^{\beta }}{\Gamma \left(\beta \right)}{\theta }^{\beta -1}{\text{e}}^{-\lambda \theta }\text{d}\lambda \text{d}\beta }}$ (17)

定理5对于对数伽玛分布(1)，若参数 $\theta$ 的先验分布密度函数由(17)式给出，则在熵损失函数下，参数 $\theta$ 的多层Bayes估计为

${\delta }_{HB}\left(x\right)=\frac{{\displaystyle {\int }_0^1{\displaystyle {\int }_0^c\frac{{\lambda }^{\beta }}{\Gamma \left(\beta \right)}\frac{\Gamma \left(n\alpha +\beta \right)}{{\left(T+\lambda \right)}^{n\alpha +\beta }}\text{d}\lambda \text{d}\beta }}}{{\displaystyle {\int }_0^1{\displaystyle {\int }_0^c\frac{{\lambda }^{\beta }}{\Gamma \left(\beta \right)}\frac{\Gamma \left(n\alpha +\beta -1\right)}{{\left(T+\lambda \right)}^{n\alpha +\beta -1}}\text{d}\lambda \text{d}\beta }}}$ (18)

证明 由θ的先验分布密度函数(17)式，所以θ的后验密度函数为

$\begin{array}{c}h\left(\theta |x\right)=\frac{\pi \left(\theta \right)f\left(x_1,\cdots ,x_n;\theta ,\alpha \right)}{{\displaystyle {\int }_0^{\infty }\pi \left(\theta \right)f\left(x_1,\cdots ,x_n;\theta ,\alpha \right)\text{d}\theta }}=\frac{{\displaystyle {\int }_0^1{\displaystyle {\int }_0^c\frac{{\lambda }^{\beta }}{\Gamma \left(\beta \right)}{\theta }^{n\alpha +\beta -1}{\text{e}}^{-\left(T+\lambda \right)\theta }\text{d}\lambda \text{d}\beta }}}{{\displaystyle {\int }_0^{\infty }{\displaystyle {\int }_0^1{\displaystyle {\int }_0^c\frac{{\lambda }^{\beta }}{\Gamma \left(\beta \right)}{\theta }^{n\alpha +\beta -1}{\text{e}}^{-\left(T+\lambda \right)\theta }\text{d}\lambda \text{d}\beta \text{d}\theta }}}}\\ =\frac{{\displaystyle {\int }_0^1{\displaystyle {\int }_0^c\frac{{\lambda }^{\beta }}{\Gamma \left(\beta \right)}{\theta }^{n\alpha +\beta -1}{\text{e}}^{-\left(T+\lambda \right)\theta }\text{d}\lambda \text{d}\beta }}}{{\displaystyle {\int }_0^1{\displaystyle {\int }_0^c\frac{{\lambda }^{\beta }}{\Gamma \left(\beta \right)}\frac{\Gamma \left(n\alpha +\beta \right)}{{\left(T+\lambda \right)}^{n\alpha +\beta }}\text{d}\lambda \text{d}\beta }}}\end{array}$ (19)

其中 $T={\displaystyle \underset{i=1}{\overset{n}{\sum }}\ln X_i}$。

$\begin{array}{c}E\left({\theta }^{-1}|X\right)=\frac{{\displaystyle {\int }_0^{\infty }{\displaystyle {\int }_0^1{\displaystyle {\int }_0^c\frac{{\lambda }^{\beta }}{\Gamma \left(\beta \right)}{\theta }^{n\alpha +\beta -2}{\text{e}}^{-\left(T+\lambda +a\right)\theta }\text{d}\lambda \text{d}\beta \text{d}\theta }}}}{{\displaystyle {\int }_0^1{\displaystyle {\int }_0^c\frac{{\lambda }^{\beta }}{\Gamma \left(\beta \right)}\frac{\Gamma \left(n\alpha +\beta \right)}{{\left(T+\lambda \right)}^{n\alpha +\beta }}\text{d}\lambda \text{d}\beta }}}\\ =\frac{{\displaystyle {\int }_0^1{\displaystyle {\int }_0^c\frac{{\lambda }^{\beta }}{\Gamma \left(\beta \right)}\frac{\Gamma \left(n\alpha +\beta -1\right)}{{\left(T+\lambda \right)}^{n\alpha +\beta -1}}\text{d}\lambda \text{d}\beta }}}{{\displaystyle {\int }_0^1{\displaystyle {\int }_0^c\frac{{\lambda }^{\beta }}{\Gamma \left(\beta \right)}\frac{\Gamma \left(n\alpha +\beta \right)}{{\left(T+\lambda \right)}^{n\alpha +\beta }}\text{d}\lambda \text{d}\beta }}}\end{array}$ (20)

3. 数值模拟

下面用R软件进行数值模拟：取形状参数 $\alpha =2$，尺度参数 $\theta =0.82$ 的对数伽玛分布 $n=50$ 个随机数，由式 $T={\displaystyle \underset{i=1}{\overset{n}{\sum }}\ln x_i}$ 计算可得 $T=79.7130$。分别根据定理1、定理2计算先验分布为Quasi先验分布时，不同损失函数下对数伽玛分布参数θ的Bayes估计值 ${\delta }_B$。

1) Linex损失函数下Quasi分布为先验分布时，尺度参数θ的Bayes估计值 ${\delta }_B$。见表1：

Table 1. Parameter estimation of parameter theta under Linex loss function

表1. Linex损失下参数θ的参数估计

2) 熵损失函数下Quasi分布为先验分布时，尺度参数θ的Bayes估计值 ${\delta }_B$。见表2：

Table 2. Parameter estimation of parameter theta under entropy loss function

表2. 熵损失下参数θ的参数估计

根据表1和表2中的极值和均方误差(MSE)，可以看出：Linex损失函数下先验分布为Quasi分布时尺度参数θ的Bayes估计值的稳健性远比熵损失函数下的尺度参数的稳健性好，横向极差不超过0.0105，纵向极差不超过0.0066。从统计决策中稳健性角度考虑，第一个条件下的 ${\delta }_B$ 更稳健。

另外，由偏差 $\Delta \delta =\left|{\delta }_B-{\theta }_0\right|$ (其中 ${\delta }_B$ 为参数θ的估计量， ${\theta }_0$ 为参数θ的真值)得到表1的偏差区间为[0.0002, 0.0081]，可见第一个条件下的偏差很小，所以精确度较高。

3) Linex损失函数下，先验分布分别为Quasi分布和伽玛分布时，尺度参数θ的Bayes估计值 ${\delta }_B$ 的对比：

由文献 [3] 可知，Linex损失函数下先验分布为伽玛分布的对数伽玛分布尺度参数θ的Bayes估计(表3)：

Table 3. Parameter estimation of parameter theta under the prior distribution of gamma distribution

表3. 伽玛分布先验分布下参数θ的参数估计

注：n = 50, T = 79.7130, α = 2。

比较上述表1和表3，表中的极值和均方误差，Linex损失函数下适当选择参数 $a,d$，对数伽玛分布尺度参数θ在Quasi先验分布下的Bayes估计比在伽玛分布为先验分布下的估计值更接近真值，更稳健。

4. 结论

数值模拟的结果表明，尺度参数θ的先验分布为Quasi分布时，Linex损失函数下的Bayes估计稳健性较熵损失函数下的稳健性更好，Bayes估计值更接近真值。另外，与文献 [3] 在同样的Linex损失函数下，将先验分布由伽玛分布改为Quasi先验分布时作对比，表1和表3中极值和均方误差结果表明适当选择Linex损失函数及Quasi先验分布的参数 $a,d$ 时，对数伽玛分布尺度参数θ在Quasi先验分布下的Bayes估计比在伽玛分布为先验分布下的估计值更接近真值，更稳健。

基金项目

国家自然科学基金项目(11801488)；新疆师范大学教学研究与改革项目(SDJG2020-30)；新疆师范大学科研发展专项项目(XJNUZX202001)。

参考文献

NOTES

^*第一作者。

^#通讯作者。