Change of variables in a Gaussian integral in matrix form

Question

I have a problem in which I have to compute the following integral: $$\mathop{\idotsint\limits_{\mathbb{R}^k}}_{\sum_{i=1}^k y_i=x} e^{-N^2r(\sum_{i=1}^k y_i^2-\frac{1}{k}x^2)} dy_1\dots dy_k,$$ where this notation means that I want to integrate over $\mathbb{R}^k$ restricted to the plane where $\sum_{i=1}^k y_i=x$ (a convolution of gaussians) and $N$ and $r$ are positive real constants. I have tried two different methods for computing this integral, but they are yielding different results. I would appreciate it very much if someone could take a look and tell me what I'm doing wrong.

Method 1

In method 1 I just wrote it as $$\mathop{\idotsint\limits_{\mathbb{R}^k}}_{\sum_{i=1}^{k}y_i=x} e^{-N^2r(\sum_{i=1}^{k}y_i^2-\frac{1}{k}x^2)} dy_1\dots dy_k =\int_{-\infty}^{\infty}\dots\int_{-\infty}^\infty e^{-N^2r((x-y_1)^2+\sum_{i=1}^{k-2}(y_i-y_{i+1})^2+y_{k-1}^2-\frac{1}{k}x^2)} \, dy_1\dots dy_{k-1}=\sqrt{\frac{1} {\pi r^{k-1}k}} \frac{\pi^k}{N^{k-1}}$$

I deduced this formula by induction, first integrating in $y_{k-1}$, then $y_{k-2}$ and so on.

Method 2

In method 2 I tried writting the function in a matrix form $$\mathop{\idotsint\limits_{\mathbb{R}^k}}_{\sum_{i=1}^{k}y_i=x} e^{-N^2r(\sum_{i=1}^{k}y_i^2-\frac{1}{k}x^2)} dy_1\dots dy_{k}=\mathop{\idotsint\limits_{\mathbb{R}^k}}_{\sum_{i=1}^{k}y_i=x} e^{-N^2r(\vec{y},Q\vec{y})} dy_1\dots dy_{k}$$ where \begin{equation} Q:=\left(\begin{array}{cccccccc} (1-\frac{1}{k})& -\frac{1}{k} & -\frac{1}{k} & \cdots & -\frac{1}{k} \\ -\frac{1}{k} & (1-\frac{1}{k}) & -\frac{1}{k} & \cdots & -\frac{1}{k} \\ \vdots & \ddots & & &\vdots \\ -\frac{1}{k} & \dots & &-\frac{1}{k} &(1-\frac{1}{k}) \end{array}\right). \end{equation}

This matrix $Q$ has eigenvalues $\lambda_0=0$, $\lambda_l=1$ and corresponding normalized eigenvetors \begin{equation} \vec{\lambda}_l=\frac{1}{\sqrt{k}}\left(\begin{array}{c} 1 \\ e^{\frac{2\pi i}{k}1l} \\ \vdots \\ e^{\frac{2\pi i}{k}(k-1)l} \end{array}\right) \end{equation} for $0\le l\le k-1$.

As I understand it, the restriction in the integral means that I shouldn't integrate in the $\lambda_0$ direction, since in this direction I must have all components equal, and the only place where the components are equal and the bound is satisfied is $(\frac{x}{k},\dots,\frac{x}{k})$. So my integration should occour in the orthogonal complement of this vector, which is a hyperplane of dimension $k-1$. Everything seems to check to this point, so I diagonalized the matrix $Q=U\Lambda U^{-1}$ and so

$$(\vec{y},Q\vec{y})=(\vec{\xi},\Lambda\vec{\xi})=\sum_{i=1}^{k-1}\xi_i^2.$$

The change of variables $\vec{\xi}=U^{-1}\vec{y}$ has a Jacobian $\frac{1}{\sqrt{k^{k-1}}}$, since $U^{-1}$ is the DFT matrix times $\frac{1}{\sqrt{k^{k-1}}}$ and the DFT matrix is known to be unitary. So

$$\mathop{\idotsint\limits_{\mathbb{R}^k}}_{\sum_{i=1}^{k}y_i=x} e^{-N^2r(\vec{y},Q\vec{y})} dy_1\dots dy_{k}=\idotsint\limits_{\mathbb{R}^k} e^{-N^2r\sum_{i=1}^{k-1}\xi_i^2} \frac{1}{\sqrt{k^{k-1}}}d\xi_1\dots d\xi_{k-1}= \sqrt{\frac{\pi^{k-1}}{k^{k-1}r^{k-1}}}\frac{1}{N^{k-1}}.$$

These two results are different and I cannot figure out why.

Thank you all in advance for your help!

Answer 1

$\newcommand\R{\mathbb R}\newcommand\1{\mathbf1}$When you say "I want to integrate over $\mathbb{R}^k$ restricted to the plane where $\sum_{i=1}^{k}y_i=x$", you have to specify the measure over the plane over which you want to integrate.

It appears you want this measure to be induced by the Lebesgue measure on $\R^k$. Then the integration can be done as follows. Let $$c:=N^2r\in(0,\infty)$$ and $$t:=x/\sqrt k,$$ the (signed) distance from the the origin to your plane $$\Pi_t:=\{y\in\R^k\colon u\cdot y=t\}=\{y\in\R^k\colon \1\cdot y=x\},$$ where $\cdot$ denotes the dot product, $\1:=(1,\dots,1)\in\R^k$, and $$u:=\1/\sqrt k$$ is a unit normal vector to the plane $\Pi_t$. Thus, instead of the parameter $x$, we use the more geometrical parameter $t$.

Then the integral in question can be written as $$I_t:=e^{ct^2}J_t,\quad\text{where}\quad J_t:=\int_{\Pi_t}\mu_t(dy)e^{-c|y|^2},$$ $|y|$ is the Euclidean norm of $y$, and, for each real $t$, $\mu_t$ is the measure over the plane $\Pi_t$ induced by the Lebesgue measure on $\R^k$ in the following sense: \begin{equation} \int_a^b dt\, \int_{\Pi_t}\mu_t(dy)g(y) =\int_{\Pi_{a,b}}dy\,g(y) \tag{1} \end{equation} for all nonnegative Borel-measurable functions $g\colon\R^k\to\R$ and all real $a$ and $b$ such that $a<b$, where $$\Pi_{a,b}:=\bigcup_{t\in[a,b]}\Pi_t=\{y\in\R^k\colon a\le u\cdot y\le b\}.$$

Then for such $a$ and $b$ we have $$\int_a^b dt\, J_t=\int_a^b dt\, \int_{\Pi_t}\mu_t(dy)e^{-c|y|^2} =K_{a,b}:=\int_{\Pi_{a,b}}dy\,e^{-c|y|^2}.$$ (See the remark on this at the end of this answer.)

To compute the integral $K_{a,b}$, let us use a substitution of the form $y=Qz$, where $Q$ is any orthogonal $k\times k$ matrix whose first column is the unit vector $u$, so that $y=Qz$ implies $z_1=u\cdot y$, where $z_j$ is the $j$'s coordinate of the vector $z$; such an orthogonal matrix $Q$ exists. Then we can write $$K_{a,b}=\int_{\R^k}dy\,e^{-c|y|^2}\,1(a\le u\cdot y\le b) \\ =\int_{\R^k}dz\,e^{-c|z|^2}\,1(a\le z_1\le b) \\ =\int_a^b dz_1\,e^{-cz_1^2}\int_{\R^{k-1}}dw\,e^{-c|w|^2} \\ =\int_a^b dz_1\,e^{-cz_1^2}\,(\pi/c)^{(k-1)/2}. $$ So, $$J_t=\frac d{dt}\,K_{a,t}=e^{-ct^2}\,(\pi/c)^{(k-1)/2}.$$ Thus, the integral in question is $$I_t:=e^{ct^2}J_t=(\pi/c)^{(k-1)/2}=\Big(\frac\pi{N^2r}\Big)^{(k-1)/2}.$$

This differs from both of your answers -- but you never defined the measure over which you integrate.

---

Remark: Intuitively, think of $a$ and $b$ as being close to a real number $t$, and hence to each other. We approximate the integral $\int_{\Pi_t}\mu_t(dy)e^{-c|y|^2}$ over the plane $\Pi_t$ by $\frac1{b-a}\,\int_{\Pi_{a,b}}dy\, e^{-c|y|^2}$, that is, by the integral over the thin layer $\Pi_{a,b}$ between two parallel planes $\Pi_a$ and $\Pi_b$ (close to the plane $\Pi_t$) divided by thickness $b-a$ of the layer.

Formally, we are dealing here with disintegration of a measure. That linked theorem deals only with probability measures, but it is trivially extended to finite measures. If we forget about this finiteness condition for a moment, then in that linked theorem we can choose $X=\R$, $Y=\R^k$, let the map $\pi\colon Y\to X$ be the projection map defined by $\pi(y):=u\cdot y$ for all $y\in Y=\R^k$, and let $\mu:=\lambda_k$ and $\nu:=\lambda_1$, where $\lambda_k$ is the Lebesgue measure over $\R^k$. I think the finiteness condition in that linked theorem is inessential, and the proof will hold for any Borel measures $\mu$, at least if $\mu$ is $\sigma$-finite. Alternatively, one can approximate here the Lebesgue measure over $\R^k$ by the finite Lebesgue measures over big cubes in $\R^k$.

Answer 2

The desired integral $I$ can be written as $^\ast$

$$I=e^{N^2rx^2/k}J(x),\;\;J(x)=\mathop{\idotsint}\delta\left(x-{\textstyle{\sum_{j=1}^{k}}y_j}\right) e^{-N^2r\sum_{j=1}^{k}y_j^2} dy_1\dots dy_{k}.$$ Fourier transform $J(x)$, $$F(q)=\int_{-\infty}^\infty e^{iqx}J(x)\,dx=\mathop{\idotsint} e^{-N^2r\sum_{j=1}^{k}y_j^2}e^{\sum_{j=1}^k qy_j} dy_1\dots dy_{k}=N^{-k}(\pi/r)^{k/2}e^{-kq^2/(4N^2r)}.$$ Now Fourier transform back from $F(q)$ to $J(x)$ and you're done: $$J(x)=(2\pi)^{-1}\int_{-\infty}^\infty e^{-iqx}F(q)\,dq=k^{-1/2} (N^2 r/\pi)^{(1-k)/2} e^{-N^2 r x^2/k}$$ $$\Rightarrow I=k^{-1/2} (N^2 r/\pi)^{(1-k)/2}.$$ This differs from both of the answers in the OP, but answer number 1 will agree if a typo is corrected ($\pi^k$ should be $\pi^{k/2}$, which is the factor coming from the Gaussian integrals in method 1).

---

$^\ast$ As noted by Iosif Pinelis, the notation $\idotsint_{\sum_{i=1}^{k}y_i=x}$ in the OP is ambiguous, I replaced it by the delta function $\idotsint \,\delta(x-\sum_{i=1}^{k}y_i)$ to remove the ambiguity.
More generally, one can define $I_a$ and $J_a(x)$ with the delta function $\delta(ax-a\sum_{i=1}^{k}y_i)$ for any $a>0$. This amounts to the same Euclidian measure on the hyperplane, but the distance from the hyperplane is rescaled by a factor $a$. Since $\delta(au)=a^{-1}\delta(u)$, one has $$I_a=a^{-1}I_1(x)=a^{-1}k^{-1/2}(N^2r/\pi)^{(1-k)/2}.$$ Iosif's answer corresponds to the choice $a=k^{-1/2}$.

Answer 3

$\newcommand\1{\mathbf1}\newcommand{\R}{\mathbb R}\newcommand{\la}{\lambda}$Here is a more explicit way to define the disintegration of the Lebesgue measure over $\R^k$ into the measures $\mu_t$ over the planes $\Pi_t$ introduced in my other answer on this page, of course with the same result.

Let us recall notations introduced in that answer: $$c:=N^2r\in(0,\infty),\quad t:=x/\sqrt k,$$ $$\Pi_t:=\{y\in\R^k\colon u\cdot y=t\}=\{y\in\R^k\colon \1\cdot y=x\},$$ where $\cdot$ denotes the dot product, $\1:=(1,\dots,1)\in\R^k$,
$$u:=\1/\sqrt k;$$

The integral in question still is \begin{equation*} I_t:=e^{ct^2}J_t,\quad\text{where}\quad J_t:=\int_{\Pi_t}\mu_t(dy)e^{-c|y|^2}, \tag{1} \end{equation*} $|y|$ is the Euclidean norm of $y$, and (for each real $t$) $\mu_t$ is the measure over the plane $\Pi_t$ now explicitly defined as follows.

Let $T\colon\R^{k-1}\to\Pi_0$ be any linear isometry of $\R^{k-1}$ onto the $(k-1)$-dimensional linear subspace $\Pi_0$ of $\R^k$. Any such isometry $T$ is given by the formula $Tz=\sum_{j=1}^{k-1}z_jb_j$ for all $z=(z_1,\dots,z_{k-1})\in\R^{k-1}$, where $(b_1,\dots,b_{j-1})$ is any orthonormal basis of $\R^{k-1}$. For each real $t$, we have $\Pi_t=\Pi_0+tu$, and so, we can define the affine isometry $U_t\colon\R^{k-1}\to\Pi_t$ of $\R^{k-1}$ onto the $(k-1)$-dimensional affine subspace $\Pi_t$ by the formula \begin{equation*} U_tz:=Tz+tu \end{equation*} for all $z\in\R^{k-1}$. Let now \begin{equation*} \mu_t:=\mu^T_t:=\la_{k-1}U_t^{-1}, \tag{2} \end{equation*} the pushforward measure for the Lebesgue measure $\la_{k-1}$ over $\R^{k-1}$ under the isometry $U_t$, so that $\mu_t(B)=\la_{k-1}(U_t^{-1}(B))=\la_{k-1}(T^{-1}(B-tu))$ for all Borel subsets $B$ of $\Pi_t$.

Remark: The measures $\mu_t=\mu^T_t$ do not depend on the choice of a linear isometry $T$ of $\R^{k-1}$ onto $\Pi_0$. Indeed, if $T$ is any such isometry, then any other such isometry is of the form $S:=TQ$, where $Q$ is any linear isometry of $\R^{k-1}$ (onto $\R^{k-1}$). Hence, \begin{equation*} \mu^S_t(B)=\la_{k-1}(S^{-1}(B-tu))=\la_{k-1}(Q^{-1}(T^{-1}(B-tu))) =\la_{k-1}(T^{-1}(B-tu))=\mu^T_t(B) \end{equation*} for all Borel subsets $B$ of $\Pi_t$; the penultimate equality in the above display holds because the Lebesgue measure is invariant with respect to linear isometries. $\Box$

It follows from (1) and (2) that \begin{equation} J_t=\int_{\R^{k-1}}dz\,\exp\{-c|U_tz|^2\}, \end{equation} by the change-of-variables formula for the pushforward measures. Since $U_tz=Tz+tu$, $Tz\perp u$, $T$ is an isometry, and $|u|=1$, we have $|U_tz|^2=|Tz|^2+t^2|u|^2=|z|^2+t^2$. So, \begin{equation} J_t=e^{-ct^2}\int_{\R^{k-1}}dz\,\exp\{-c|z|^2\} =e^{-ct^2}(\pi/c)^{(k-1)/2}. \end{equation} Thus, the integral in question is $$I_t=e^{ct^2}J_t=(\pi/c)^{(k-1)/2}=\Big(\frac\pi{N^2r}\Big)^{(k-1)/2},$$ which is what was obtained a bit differently in my other answer.

Answer 4

$\newcommand\1{\mathbf1}\newcommand{\R}{\mathbb R}\newcommand{\la}{\lambda}$Here is yet another, "multivariate calculus" way to treat your integral, with the same result as in my other two answers on this page.

Let $$c:=N^2r\in(0,\infty),$$ $$\Pi_x:=\{y\in\R^k\colon \1_k\cdot y=x\},$$ where $\cdot$ denotes the dot product and $\1_k:=(1,\dots,1)\in\R^k$.

The integral in question is \begin{equation*} I_x:=e^{cx^2/k}J_x,\quad\text{where}\quad J_x:=\int_{\Pi_x}\mu_x(dy)e^{-c|y|_k^2}, \tag{1} \end{equation*} $|y|_k$ is the Euclidean norm of $y\in\R^k$, and $\mu_x(dy)$ is the surface area element on the plane $\Pi_x$. This plane is the graph of the function $f_x\colon\R^{k-1}\to\R$ given by the formula $f_x(y):=x-y_1-\cdots-y_{k-1}$ for $y=(y_1,\dots,y_{k-1})\in\R^{k-1}$. So (cf. the case $k=3$), \begin{align*} J_x&=\int_{\R^{k-1}}dy_1\cdots dy_{k-1}\,\sqrt{1+\sum_{j=1}^{k-1}\Big(\frac{\partial f_x}{\partial y_j}\Big)^2} \\ & \times \exp\Big\{-c(x-y_1-\cdots-y_{k-1})^2-c\sum_{j=1}^{k-1}y_j^2\Big\} \\ &=\sqrt k\,\int_{\R^{k-1}}dy\, \exp\big\{-c(x-y\cdot\1_{k-1})^2-c|y|_{k-1}^2\big\} \\ &=\sqrt k\,\int_{\R^{k-1}}dz\, \exp\big\{-c(x-z_1\sqrt{k-1})^2-c|z|_{k-1}^2\big\}; \end{align*} for the last displayed equality, we use any substitution of the form $y=Qz$, where $Q$ is any orthogonal $(k-1)\times(k-1)$ matrix whose first column is $\1_{k-1}/\sqrt{k-1}$ and $z=(z_1,\dots,z_{k-1})$, so that $z_1=y\cdot\1_{k-1}/\sqrt{k-1}$.

So, \begin{align*} J_x&=\sqrt k\,\int_\R dz_1\, \exp\big\{-c(x-z_1\sqrt{k-1})^2-cz_1^2\big\} \\ &\times\int_{\R^{k-2}}dz_2\cdots dz_{k-1}\, \exp\Big\{-c\sum_{j=2}^{k-1}z_j^2\Big\} \\ &=\sqrt k\ \;\frac{e^{-cx^2/k}(\pi/c)^{1/2}}{\sqrt k}\ \;(\pi/c)^{(k-2)/2} \\ &=e^{-cx^2/k}\;(\pi/c)^{(k-1)/2}. \end{align*}

Thus, by (1), the integral in question is $$I_x=(\pi/c)^{(k-1)/2}=\Big(\frac\pi{N^2r}\Big)^{(k-1)/2},$$ which is what was obtained a bit differently in my other two answers.