Any exponential family distribution can be expressed as $p(x|\theta) = g(\theta) f(x) e^{\phi(\theta)^T T(x)} = f(x) e^{\phi(\theta)^T T(x) - A(\theta)}$ where $A(\theta) = -\log{g(\theta)}$.

We know that this integrates to 1 since it's a valid pdf:

$\int f(x) e^{\phi(\theta)^T T(x) - A(\theta)} dx = 1$

The expected value of the sufficient statistic is obtained by differentiating both sides with respect to $\theta$. Assuming this is a parameter vector, we will need a vector derivative:

$\begin{align} \frac{\partial}{\partial \theta_i} \int f(x) e^{\phi(\theta)^T T(x) - A(\theta)} dx &= \frac{\partial}{\partial \theta_i} (1) \\ \int \frac{\partial}{\partial \theta_i} [ f(x) e^{\phi(\theta)^T T(x) - A(\theta)}] dx &= 0 \\ \int p(x|\theta) \left[ \frac{\partial \phi(\theta)^T}{\partial \theta_i} T(x) - \frac{\partial A(\theta)}{\partial \theta_i } \right] dx &= 0 \\ \nabla_i \phi(\theta)^T \int p(x|\theta) T(x) dx - \nabla_i A(\theta) \int p(x|\theta) dx &= 0 \\ \nabla_i \phi(\theta)^T E \left[ T(x) \right] - \nabla_i A(\theta) &= 0\end{align}$

But this is not quite correct since the first term is a scalar (dot product between two vectors?) and the second term is a vector.

The expectation should be the vector i.e. I should have something like $E[T(x)]_i$ with the vector index $i$ on the expected value. What has gone wrong?

I think (but not sure) that the answer should be $E[T(x)]_i = \frac{1}{\phi'(\theta)} \nabla_iA(\theta)$


1 Answer 1


Let $\theta \in \Theta \subset \mathbb R^m$ and $\phi : \Theta \to \mathbb R^d$ so $T(x) \in \mathbb R^d$ too. I'll use without proof the fact that differentiation and integration can be exchanged for exponential families.

First I'm going to treat this as a natural exponential family and see what happens (this is like if $\phi$ is just the identity function). I'll use $\theta$ as the parameter so right now I have $p(x|\theta) = \exp(\theta^TT(x) - A(\theta))f(x)$.

In this case the differentiation is easy and we have $$ 0 = \int p(x|\theta) \left[T_j(x) - \nabla _jA(\theta)\right]\,\text dx $$ so $$ \nabla A(\theta) = \text E(T(X)). $$

So if you're content to work with $\phi := \phi(\theta)$ as the parameter, then this is a tidy result.

Now I'll consider what happens if $m = d = 1$ to build our intuition for the general case. We'll have $$ 0 = \int p(x|\theta)\left[T(x) \cdot \phi'(\theta) - A'(\theta) \right]\,\text dx \\ \implies \text E(T(X)) \cdot \phi'(\theta) = A'(\theta) $$ so $$ \text E(T(X)) = \frac{A'(\theta)}{\phi'(\theta)} $$ as you thought.

But the situation is more complicated if $m,d > 1$. The derivates in particular are more complicated because now $\phi$ is vector-valued so we'll get a matrix of first derivatives.

We need $$ \frac{\partial}{\partial \theta_j} \phi(\theta)^TT(x). $$ Writing this as a sum, we have $$ \frac{\partial}{\partial \theta_j} \phi(\theta)^TT(x) = \sum_{i=1}^d T_i(x) \frac{\partial}{\partial \theta_j}\phi_i(\theta). $$ $\phi : \Theta\to\mathbb R^d$ with $\Theta\subset\mathbb R^m$ so for each component we have $\phi_i : \Theta\to\mathbb R$ with $i=1,\dots,d$. Each $\phi_i$ is just a typical scalar-valued function so $\frac{\partial}{\partial \theta_j} \phi_i = (\nabla \phi_i)_j$.

I'll collect these gradients into a matrix $\Phi(\theta) \in \mathbb R^{m\times d}$ where $$ \Phi(\theta) = \left[\begin{array}{c|c|c|c}\nabla \phi_1(\theta) & \nabla \phi_2(\theta) & \cdots & \nabla \phi_d(\theta)\end{array} \right]. $$

This means $$ \frac{\partial}{\partial \theta_j}\phi_i(\theta) = \Phi(\theta)_{ji} $$ so $$ \frac{\partial}{\partial \theta_j} \phi(\theta)^TT(x) = \left(\Phi(\theta)T(x)\right)_j. $$

All together this means $$ \mathbf 0 = \int p(x|\theta) \left[ \Phi(\theta)T(x) - \nabla A(\theta) \right] \,\text dx \\ \implies \Phi(\theta)\text E [T(X)] = \nabla A(\theta) $$ which is a linear system that we're trying to solve for $\text E [T(X)]$. The dimensions in question and the specific properties of $\phi$ now will determine how many solutions this system has.

Here's one example. Suppose $X \sim \mathcal N(\mu, \sigma^2)$ so $$ p(x|\mu,\sigma^2) = \frac{1}{\sqrt{2\pi\sigma^2}}\exp\left(-\frac{(x-\mu)^2}{2\sigma^2}\right). $$ Rearranging and expanding the quadratic we have $$ p(x|\mu,\sigma^2) = \frac{1}{\sqrt{2\pi}}\exp\left(\frac{\mu}{\sigma^2}x - \frac{1}{2\sigma^2} x^2 - \left[\frac{\mu^2}{2\sigma^2} + \frac 12 \log \sigma^2\right]\right) $$ which means we'll have $\theta = (\mu,\sigma^2)$, $\Theta = \mathbb R \times (0,\infty)$, $T(x) = {x \choose x^2}$, $$ \phi(\theta) = \left(\begin{array}{c}\frac{\theta_1}{\theta_2} \\ \frac{-1}{2\theta_2}\end{array}\right), $$ and $$ A(\theta) = \frac 12 \theta_1^2\theta_2^{-1} + \frac 12 \log\theta_2. $$

This means $$ \Phi(\theta)= \left(\begin{array}{cc}\theta_2^{-1} & 0 \\ -\theta_1\theta_2^{-2} & \frac 12 \theta_2^{-2}\end{array}\right). $$ Noting that $$ \det \Phi(\theta) = \frac 12 \theta_2^{-3} > 0 $$ we know this matrix is invertible, and it's just $2\times 2$ so the inverse is easy enough to compute: $$ \Phi^{-1}(\theta) = 2\theta_2^3 \left(\begin{array}{cc}\frac 12 \theta_2^{-2} & 0 \\ \theta_1\theta_2^{-2} & \theta_2^{-1}\end{array}\right) = \left(\begin{array}{cc}\theta_2 & 0 \\ 2\theta_1\theta_2 & 2\theta_2^2\end{array}\right). $$ We also have $$ \nabla A(\theta) = \left(\begin{array}{cc}\theta_1\theta_2^{-1} \\ -\frac 12\theta_1^2\theta_2^{-2} + \frac 12 \theta_2^{-1}\end{array}\right) $$ so $$ \text E(T(X)) = \Phi^{-1}(\theta)\nabla A(\theta) \\ = \left(\begin{array}{c}\theta_1 \\ \theta_1^2 + \theta_2\end{array}\right) = {\mu \choose \mu^2 + \sigma^2} $$ which we know to be the first two raw moments of $X$. This was a lot of work to get there but it shows the unavoidable role of the matrix $\Phi$.

Update: here's the multivariate normal case. Let $X \sim \mathcal N_n(\mu, \Omega^{-1})$ where $\Omega$ is the precision matrix. I'll use $\Omega_{(i)}$ for the $i$th row of $\Omega$, and for a vector $v$ and $1\leq i\leq j$ I'll use $v_{i:j}$ to be $(v_i,v_{i+1},\dots,v_j)$.

The density is $$ p(x \mid \mu, \Omega) = \exp\left[-\frac 12 x^T\Omega x + x^T \Omega \mu - \frac 12 \left(\mu^T\Omega\mu - \log\det \Omega\right)\right] \cdot(2\pi)^{-n/2}. $$

The parameters are $\mu$ and $\Omega$ but I want to deal with these as a flat vector so I'll combine $\mu$ with $\Omega$ "unrolled" as $$ \theta = (\mu, \Omega_{11}, \Omega_{12}, \dots, \Omega_{1n}, \Omega_{22}, \dots, \dots, \Omega_{n-1,n}, \Omega_{nn}). $$ I'm only taking the diagonal and one triangle of $\Omega$ since it's symmetric. This means $\Theta\subset \mathbb R^m$ where $m = n + {n+1\choose 2} = 2n + {n\choose 2}$.

I now need to get this as $\phi(\theta)^TT(x) - A(\theta)$. For $\phi$ and $T$ I'll use the fact that $$ x^T\Omega x = 2\sum_{i < j} x_ix_j \Omega_{ij} + \sum_{i=1}^n x_i^2\Omega_{ii} $$ so I'll take $$ \phi(\theta) = \left(\Omega\mu, -\frac 12 \Omega_{11}, -\Omega_{12}, \dots, -\Omega_{1n}, -\frac 12 \Omega_{22}, -\Omega_{23}, \dots, \dots, -\Omega_{n-1,n}, -\frac 12 \Omega_{nn}\right) $$ (I'm writing this in terms of $\mu$ and $\Omega$ instead of $\theta$ just because it's way clearer, but these can all be written in terms of $\theta$ if desired) and $$ T(x) = (x, x_1^2, x_1x_2, \dots, x_1x_n, x_2^2, x_2x_3, \dots,\dots, x_n^2). $$

We can now construct $\Phi$. I'll build it as a $2\times 2$ block matrix where $$ \Phi(\theta) = \left(\begin{array}{c|c}A&B\\\hline C&D\end{array}\right) $$ with $A$ being $n\times n$ and $D$ being ${n+1c\choose 2}\times{n+1\choose 2}$. $A$, $B$, and $D$ are easy to work out.

For $A$, I need $\frac{\partial \phi_i}{\partial \theta_j}$ where $i,j=1\dots,n$. Since $\phi_i(\theta) = \Omega_{(i)}^T\mu$ for this range of $i$ values, and $\theta_j = \mu_j$, we'll get $$ (\nabla \phi_i)_{1:n} = \Omega_{(i)} $$ (using the symmetry of $\Omega$ for the fact that $\Omega_{k\ell} = \Omega_{\ell k}$). This means that $A = \Omega$.

Now for $i > n$, $\phi_i(\theta)$ is proportional to some $\Omega_{j\ell}$ so derivatives w.r.t. $\mu$ are zero. This means $B = \mathbf O$, the appropriately sized zero matrix.

Next, for $i,j>n$ if $i\neq j$ then $\frac{\partial \phi_i}{\partial \theta_j} = 0$ so $D$ is diagonal. If $i$ corresponds to a diagonal entry of $\Omega$ then we have a derivative of $-\frac 12$, otherwise it'll be a derivative of $-1$. The first element of the diagonal of $D$ corresponds to $\Omega_{11}$ so it begins with $-\frac 12$. The next value of $-\frac 12$ is for $\Omega_{22}$ which happens at element $1+n$ since all of the first row of $\Omega$ has been counted and $\Omega_{21}$ is skipped (since it's accounted for by $\Omega_{12}$). There are $n-2$ remaining elements of $\Omega_{(2)}$ to account for, so the 3rd value of $-\frac 12$ in $D$ is at position $1+n+(n-1)$. The next $-\frac 12$ is at $1+n+(n-1)+(n-2)$ and so on, so the locations of $-\frac12$s are given by the sequence $$ 1+\sum_{j=0}^{i-2}(n-j) $$ with $i=1,\dots,n$. As an example, if $n=5$ then the $-\frac 12$s will be at locations $1,6,10,13,15$. As a sanity check, $D$ is supposed to have a diagonal of length ${n+1\choose 2}$, and the final $-\frac 12$ is for $\Omega_{nn}$. We have $$ 1+\sum_{j=0}^{n-2}(n-j) = \sum_{j=0}^{n-1}(n-j) \\ = \sum_{j=1}^n j = {n+1\choose 2} $$ so this agrees with that.

Finally, we need $C$. I'm going to skip $C$ for now until the final step.

For $\Phi^{-1}$, $\Omega$ and $D$ are always invertible so $\Phi$ is too, and we can use the formula for the inverse of a $2\times 2$ block matrix to find $$ \Phi^{-1}(\theta) = \left(\begin{array}{c|c}\Omega^{-1} & \mathbf O \\ \hline -D^{-1}C\Omega^{-1} & D^{-1}\end{array}\right). $$

Now we can get $\nabla A$ to finish this off. I'll also partition it into the first $n$ elements, which lines up the upper two blocks of $\Phi^{-1}$, and the remaining ${n+1\choose 2}$ which I'll call $V$. We have $$ A(\theta) = \frac 12 \mu^T\Omega\mu - \frac 12 \log\det \Omega $$ so $$ \frac{\partial A}{\partial \theta_{1:n}} = \frac{\partial A}{\partial \mu} = \Omega\mu. $$ This means $$ \nabla A(\theta) = {\Omega\mu \choose V} $$ so $$ \Phi^{-1}(\theta)\nabla A(\theta) = {\mu \choose -D^{-1}(C\mu - V)}. $$ This shows that $\text E(T_{1:n}(X)) = \text E(X) = \mu$ which is nice to see.

For the rest, I won't work it out in general but I'll just show what we get for the elements corresponding to diagonal terms of $\Omega$ and off-diagonal terms.

For $V$, I'll need derivatives of $\log\det \Omega$ w.r.t. $\Omega_{ij}$ which I can get via this result: $$ \frac{\partial}{\partial \Omega_{ij}} \log\det \Omega = (2 \Sigma - \text{diag}(\Sigma)) $$ due to symmetry. This means that for the elements of $V$ corresponding to a diagonal element of $\Omega$ we'll have $\frac 12 \mu_i^2 - \frac 12 \Sigma_{ii}$ while the off-diagonal will be $\mu_i\mu_j - \Sigma_{ij}$.

For $C\mu$, the diagonal elements contribute $\mu_i^2$ because for that row the only non-zero element is $\mu_i$ in the $i$th column. $-D^{-1}$ will multiply the result by $2$, and then $V_i$ is $$ 2\left(\mu_i^2 - \frac 12 \mu_i^2 + \frac 12 \Sigma_{ii} \right) \\ = \mu_i^2 + \Sigma_{ii}. $$ For the off-diagonal, you can show that $C\mu$ will contribute $2\mu_i\mu_j$ so we get $$ 2\mu_i\mu_j - \mu_i\mu_j + \Sigma_{ij} = \mu_i\mu_j + \Sigma_{ij}. $$

Thus after all this work we have confirmed that $\text E(X_i) = \mu_i$ and $\text E(X_iX_j) = \mu_i\mu_j + \Sigma_{ij}$.

  • $\begingroup$ Wow. Thanks for an awesome reply! But what about if I want to consider e.g. a multivariate normal distribution in D dimensions. In this case the matrix of derivatives $\Phi(\theta) $ will be DxD meaning I won't be able to invert it explicitly. And yet we should still be able to obtain algebraic solution for E[T(X)] in this same way? $\endgroup$
    – user11128
    Sep 7, 2019 at 22:28
  • $\begingroup$ @user11128 I think in general $\Phi(\theta)^{-1}\nabla A(\theta)$ might be as explicit as you can get if you take this approach unless $\Phi$ happens to have some special structure. I'll go through the multivariate normal example and update with that. I'm also looking more into general conditions for $\Phi$ to be invertible. $\endgroup$
    – jld
    Sep 8, 2019 at 1:01
  • $\begingroup$ @user11128 I just updated with a multivariate normal example, $\Phi$ had enough structure that the inverse was easy to do $\endgroup$
    – jld
    Sep 12, 2019 at 16:48

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.