圆域上的二元二次多项式回归模型的D-最优设计

D-Optimal Design for Duality Quadratic Polynomial Regression Models in Circle Region

摘要: 对圆域上的二元二次多项式回归模型在遵循设计点“对称 + 均匀”的前提下，给出并证明其特定的内接正方形的顶点和圆心组成的饱和设计是D-最优设计，并给出了相应设计的最小二乘估计，是把最优设计理论应用到多元多项式回归模型的又一有益的尝试。

Abstract: For duality quadratic polynomial regression models in circle region, D-op designs were given and proved according to the criterion “Symmetry + Uniform” with four vertexes of special square and central angle in Circle, and the least squares estimates were given, It is a very useful attempt to apply optimal design theory to polynomial regression models.

1. 引言

在试验设计中，令R表示设计利益区域，X是区域R内任一点，回归模型一般形为

$y=f^{\text{T}}\left(X\right)\beta +\epsilon$

其中 $f^{\text{T}}\left(X\right)=\left(f_1\left(X\right),f_2\left(X\right),\cdots ,f_p\left(X\right)\right)$ 是由模型决定的函数向量，y是响应观测值,而 ${\beta }^{\text{T}}=\left({\beta }_1,{\beta }_2,\cdots ,{\beta }_p\right)$ 是模型中的待估参数， $\epsilon$ 是误差，通常假设 $E\left(\epsilon \right)=0$ ， $D\left(\epsilon \right)={\sigma }^2$ ， $\sigma$ 是已知的。

如果用 $M\left(\xi \right)$ 表示测度设计 $\xi$ 的信息矩阵，所谓的D-最优准则就是使得 $M\left(\xi \right)$ 的行列式达到最大，而且测度设计 ${\xi }^{\ast }$ 是D-最优设计的充分必要条件 [1] [2] 是方差函数 $d=f^{\text{T}}\left(X\right)M^{-1}\left({\xi }^{\ast }\right)f\left(X\right)\le p$ (模型中待估参数的个数)。

2. 二元二次多项式回归模型的最优设计

多元多项式回归模型可以用来处理一大类非线性问题，在应用数理统计学中占有重要地位。这里讨论的二元二次多项式回归模型

$y={\beta }_0+{\beta }_{1}x+{\beta }_{2}y+{\beta }_{11}{x}^2+{\beta }_{22}{y}^2+\epsilon$ (1)

这里因子空间是圆域 $\Re =\left\{\left(x,y\right)|x^2+y^2\le 1\right\}$ ，在模型中含有五个待估参数，我们尝试用圆域上的五个设计点去进行饱和设计，同时使得这些设计点的结构既要对称 [3] ，同时分布又很均匀 [4] ，即所谓的“对称 + 均匀”原则。

对模型(1)的饱和设计 $\xi$ 的基本思想是：先取圆周的任意内接正方形的四个顶点和圆心组成的五点设计，不妨取圆心为坐标原点，两个互相垂直的直径分别为x轴，y轴，设五个设计点分别是 $\left(a,b\right)$ ， $\left(-a,-b\right)$ ， $\left(b,-a\right)$ ， $\left(-b,a\right)$ ， $\left(0,0\right)$ 。

而且采用测度设计，记每个顶点的测度为 $\alpha$ ，圆心的测度为 $\beta$ ，即测度矩阵 $D\left(\xi \right)=diag\left(\alpha ,\alpha ,\alpha ,\alpha ,\beta \right)$ ，其中 $4\alpha +\beta =1,\alpha >0,\beta >0$ 。

可将模型(1)的函数向量改写为

$f^{\text{T}}\left(X\right)=\left(1,x^2,y^2,x,y\right)$

设饱和设计 $\xi$ 的结构矩阵为X，

$X=\left[\begin{array}{ccccc}1& a^2& b^2& a& b\\ 1& a^2& b^2& -a& -b\\ 1& b^2& a^2& b& -a\\ 1& b^2& a^2& -b& a\\ 1& 0& 0& 0& 0\end{array}\right]$ ，

相应设计 $\xi$ 的信息矩阵为 $M=X^{\text{T}}D\left(\xi \right)X$ ，利用 $a^2+b^2=1$ ，则

$M=\left(\begin{array}{ccccc}1& 2\alpha & 2\alpha & 0& 0\\ 2\alpha & 2\left(a^4+b^4\right)\alpha & 4a^2{b}^2\alpha & 0& 0\\ 2\alpha & 4a^2{b}^2\alpha & 2\left(a^4+b^4\right)\alpha & 0& 0\\ 0& 0& 0& 2\alpha & 0\\ 0& 0& 0& 0& 2\alpha \end{array}\right)$

要寻求模型的D-最优设计，即使得行列式 $\left|M\right|$ 取最大值的设计，而

$\left|M\right|=16{\alpha }^4\left(1-4\alpha \right){\left(2a^2-1\right)}^2$

令 $\{\begin{array}{l}\frac{\partial \left|M\right|}{\partial \alpha }=0\\ \frac{\partial \left|M\right|}{\partial a}=0\end{array}$ ，则有效驻点是 $\{\begin{array}{l}\alpha =\frac{1}{5}\\ a=0\end{array}$ ， $\{\begin{array}{l}\alpha =\frac{1}{5}\\ a=\pm \frac{\sqrt{2}}{2}\end{array}$ ，

这两组驻点能使得 $\left|M\right|$ 取最大值的，只有在驻点 $\alpha =\frac{1}{5}$ ， $a=0$ 时达到，此时 ${\left|M\right|}_{\max }=\frac{16}{5^5}$ ，此时 $\alpha =\beta =\frac{1}{5}$ ，说明测度设计是均匀的，设计点呈对称结构，它们是 $\left(0,1\right)$ ， $\left(0,-1\right)$ ， $\left(1,0\right)$ ， $\left(-1,0\right)$ ， $\left(0,0\right)$ 。

下面证明上述五点设计对模型(1)是D-最优设计，由前文知，只要能验证

$d=f^{\text{T}}\left(X\right)M^{-1}\left({\xi }^{\ast }\right)f\left(X\right)\le 5$ 。

此时

$M^{-1}\left({\xi }^{\ast }\right)=\left(\begin{array}{ccccc}5& -5& -5& 0& 0\\ -5& \frac{15}{2}& 5& 0& 0\\ -5& 5& \frac{15}{2}& 0& 0\\ 0& 0& 0& \frac{5}{2}& 0\\ 0& 0& 0& 0& \frac{5}{2}\end{array}\right)$ ，

$d=5-\frac{15}{2}\left(x^2+y^2\right)+\frac{15}{2}\left(x^4+y^4\right)+10x^2{y}^2$

在区域 $\Re =\left\{\left(x,y\right)|x^2+y^2\le 1\right\}$ 的极大值情况为

① 在 $\Re$ 的内部 $x^2+y^2<1$ ，由 $\{\begin{array}{l}\frac{\partial d}{\partial x}=30x^3+20x{y}^2-15x=0\\ \frac{\partial d}{\partial y}=30y^3+20x^{2}y-15y=0\end{array}$ ，则

$\{\begin{array}{l}x=0\\ y=0\\ d=5\end{array}$ ， $\{\begin{array}{l}x=0\\ y^2=\frac{1}{2}\\ d=\frac{25}{8}\end{array}$ ， $\{\begin{array}{l}{x}^2=\frac{1}{2}\\ y=0\\ d=\frac{25}{8}\end{array}$ ， $\{\begin{array}{l}{x}^2=\frac{3}{10}\\ y^2=\frac{3}{10}\\ d=\frac{11}{4}\end{array}$

② 在 $\Re$ 的边界 $x^2+y^2=1$ ，有 $d=5-5x^2\left(1-x^2\right)\le 5$ ，当且仅当在 $x=0,\pm 1$ 的时候取最大值。

由①②知道方差函数，即正定二次型d的最大值是5 (等于模型中的待估参数的个数)，而且当且仅当在上述设计点处才取得最大值，这就证明了所采用的饱和均匀的等测度设计是D最优设计。

3. 模型参数的最小二乘估计

由于模型(1)的待估参数向量为 $\stackrel{\to }{\beta }={\left({\beta }_0,{\beta }_{11},{\beta }_{22},{\beta }_1,{\beta }_2\right)}^{\text{T}}$ ，模型(1)可以记为

$y=f^{\text{T}}\left(X\right)\stackrel{\to }{\beta }+\epsilon$

由最小二乘估计公式 [5]

$\widehat{\stackrel{\to }{\beta }}={\left(X^{\text{T}}X\right)}^{-1}{X}^{\text{T}}Y=\frac{1}{5}{M}^{-1}{X}^{\text{T}}Y$ (2)

这里列向量 $Y={\left(y_1,y_2,y_3,y_4,y_5\right)}^{\text{T}}$ ，其中 $y_i\left(i=1,2,3,4,5\right)$ 表示第i次试验的响应观测值。把上述的 $M^{-1},X^{\text{T}},Y$ 的结果代入式(2)，则得到 $\stackrel{\to }{\beta }$ 的最小二乘估计为

$\widehat{\stackrel{\to }{\beta }}=\left(\begin{array}{c}{\hat{\beta }}_0\\ {\hat{\beta }}_{11}\\ {\hat{\beta }}_{22}\\ \begin{array}{c}{\widehat{\stackrel{\to }{\beta }}}_1\\ {\widehat{\stackrel{\to }{\beta }}}_2\end{array}\end{array}\right)=\left(\begin{array}{c}{y}_5\\ \frac{1}{2}\left(y_1+y_2\right)-y_5\\ \frac{1}{2}\left(y_3+y_4\right)-y_5\\ \begin{array}{l}\frac{1}{2}\left(y_1-y_2\right)\\ \frac{1}{2}\left(y_3-y_4\right)\end{array}\end{array}\right)$ 。