鲁棒实验设计（ED-最优设计）

写在前面

在中文社区，对于最优设计这块的讨论几近于无。写这篇博客也算是抛砖引玉。

论文：《Robust Experiment Design via Stochastic Approximation 》（1984）
作者：LUC PRONZATO AND ERIC WALTER

介绍了一种统计实验设计方案，利用Fisher矩阵的行列式的期望评判参数估计精度。然而，Fisher矩阵的计算需要知道未知参数的值，这实际上是不能获得的，因此需要一个标称值(nominal value of the parameter.)

ED-Optimal design

D-optimal design

在介绍ED-optimal之前，需要对D-optimal进行介绍。
最大似然估计写为：
$${\theta }_{ml}=argma{x}_{\theta }[p(y|\theta )]$$
最大似然估计是渐进正态的
$$\mathcal{N}(\theta ,M^{-1}(\theta ))$$
其中，$M(\theta )$是Fisher信息矩阵，定义为
$$M(\theta )=E_{y|\theta }\frac{\partial lnp(y|\theta )}{\partial \theta }\frac{\partial lnp(y|\theta )}{\partial {\theta }^T}dy$$
也等于$$M(\theta )=E_{y|\theta }[\frac{{\partial }^{2}lnp(y|\theta )}{\partial \theta \partial {\theta }^T}]$$Fisher矩阵针对$y$求取了条件分布的期望，因此实际上是关于$\theta$的函数。
在多元线性回归（高斯噪声）的情形下，Fisher信息矩阵的表达式是
$$M(\theta )=X^T{\Sigma }^{-1}X$$其中，X为数据矩阵（设计矩阵），或者说导数矩阵（对参数向量求导）。
可以发现，这个表达式与$\theta$没有关系了。这个表达式的推导，我会另外写文章详细探讨。
然而，在非线性的情况下，上面的表达式仍然与$\theta$有关系。其原因在于X是导数矩阵，线性情形下，任何点的导数都一致。然而，对于非线性问题，我们实际上是把模型先线性化，得到这个导数矩阵。而这个导数矩阵是在哪个点求得的呢？
实际上，应该从参数的真实值处展开，然而真值是不知道的。D-optimal的做法是，寻找一个标称值${\theta }_0$，在该值处线性化正向模型，再利用线性回归的Fisher矩阵表达式。这个${\theta }_0$则只能依靠先验知识或专家决策。

D-optimal就表达为：
$$e=argma{x}_e(det(M({\theta }_0)))$$
$e$表示实验配置。说白了就是找到一组设计变量，使得估计的Fisher矩阵的行列式最大。

ED-optimal

ED-optimal所作出的改进是考虑参数$\theta$的不确定性，实际上也就是先验分布。我们在$\theta$所有可能的值上求取期望，从而获得相比单取一个${\theta }_0$更加鲁棒的结果。ED-optimal的表达式是：
$$e^{\ast }=argma{x}_e{\int }_{\theta }det(M(\theta ))p(\theta )d\theta$$其中，$p(\theta )$我个人理解是先验分布。
另一方面，从决策理论的角度看，这实际上是某种贝叶斯风险。风险函数是负的Fisher信息的行列式。
值得注意的是，上式要区别于下式：
$$e^{\ast }=argma{x}_e\det (M({\int }_{\theta }\theta p(\theta )d\theta ))$$即先求期望再计算Fisher信息。

一个例子

令$$y=e^{-\theta x}(x\ge 0)+ϵ$$
噪声服从$ϵ\sim N(0,{\sigma }^2)$
先求出Fisher信息的表达式：
$$\begin{array}{rl}\frac{{\partial }^{2}y}{\partial {\theta }^2}& =\frac{x^2}{{\sigma }^2}y{e}^{-\theta x}-2\frac{x^2}{{\sigma }^2}{e}^{-2\theta x}\\ M(\theta ,x)& =-E_{y|\theta }[\frac{{\partial }^{2}y}{\partial {\theta }^2}]\\ & =-\frac{x^2}{{\sigma }^2}{e}^{-\theta x}{E}_{y|\theta }[y]+2\frac{x^2}{{\sigma }^2}{e}^{-2\theta x}\\ & =-\frac{x^2}{{\sigma }^2}{e}^{-2\theta x}+2\frac{x^2}{{\sigma }^2}{e}^{-2\theta x}\\ & =\frac{x^2}{{\sigma }^2}{e}^{-2\theta x}\end{array}$$注意到和参数$\theta$和设计变量$x$都有关系。
ED-optimal的表达式就是
$${\int }_{\theta }\frac{x^2}{{\sigma }^2}{e}^{-2\theta x}p(\theta )d\theta$$设置$\theta$为$[a,b]$上的均匀分布，则上面的积分为：
$$-\frac{x}{2{\sigma }^2(b-a)}(e^{-2xb}-e^{-2xa})$$
最优化的过程省略。