# Sum of squared normals: $\sum_{k=1}^{n}U_k^2 = \frac{1}{n}\sum_{k=1}^{n}a_k^2 \cdot \Gamma\left(\frac{n}{2},\frac{1}{2}\right)$

Assume $$U_k \sim \mathcal{N}(0,a_k^2)$$, where $$a_k \rightarrow c > 0$$ as $$k \rightarrow \infty$$. It follows that $$U_k^2 \sim \Gamma(\frac{1}{2}, \frac{1}{2a_k^2})$$.

I'm interested in the exact and the limiting distribution ($$n \rightarrow \infty$$) of the following sum $$S_n$$, which I conjecture to be $$S_n = \sum_{k=1}^{n}U_k^2 \stackrel{d}{\approx} \frac{1}{n}\sum_{k=1}^{n}a_k^2 \cdot \Gamma\left(\frac{n}{2},\frac{1}{2}\right)$$

I know the summation of gamma random variables can be tricky with varying scale parameters (example 1, example 2), however, the previous conjecture seems to work numerically for various $$c$$ and $$n$$ values. Maybe the derivation gets simplified due to shape being $$1/2$$ for all $$k$$? Unless it works only with my specific sequence $$a_k$$.

My question: could we arrive analytically to this result? Maybe it is a known result?

Edit: example R code

require(magrittr)
ak <- function(k){
lapply(1:k, function(z){
1/z^2 - 0.5^z/z^2
}) %>% do.call(c,.) %>% sum
}

lapply(1:100,function(k){
Uk <- rnorm(n = 1000, mean = 0, sd = ak(k))
Uk^2
}) %>% do.call(rbind,.) -> res

# sum of square, dim (1000 x 1)
apply(res,2,sum) -> r0

# Comparing the distributions:
S1 <- rgamma(1000, shape = 100/2, rate = 1/2) * mean(sapply(1:100,ak)^2)
ecdf(r0) %>% plot
ecdf(S1) %>% lines(.,col=2)

• Your "limiting distribution" either doesn't exist or is a constant, depending on $(a_k).$ You can apply the Central Limit theorem directly to the $U_k$ to obtain better information about the asymptotic distribution.
– whuber
Feb 24, 2021 at 15:10
• At stats.stackexchange.com/a/72486/919 I show how to find the exact expression. Maybe that answers your question? My answer at stats.stackexchange.com/a/350351/919 is sufficiently general that it directly applies here, too.
– whuber
Feb 24, 2021 at 15:36
• I forgot that the shape parameters all have to be integral for that partial fractions approach to work. In your case it could be used to approximate the distribution (by considering terms in disjoint pairs having nearly-equal scale parameters) but not to obtain the exact distribution.
– whuber
Feb 24, 2021 at 17:40

You might find Moment Generating Functions useful here. One important property is that the MGF of $$S_n = \sum_{i=1}^nb_iX_i$$ is the product of rescaled MGFs: $$M_{S_n}(t) = \prod_{i=1}^nM_{X_i}(b_it)$$ You can incorporate the $$1/n$$ factor and $$a_k^2$$ into the $$b_i$$. Use the fact that the MGF of the $$\Gamma(\alpha, \beta)$$ distribution is $$(1 - t/\beta)^{-\alpha}$$ defined for $$t<\beta$$. This function is continuous in $$t$$ so computing its limit as $$n\to\infty$$ should be possible. You want to compute the MGF limit because if the MGFs converge then the random variables do.