How to take the variance of a second order expansion? $\text{Var}\left[aX+bY+cXY+mX^2+nY^2\right]$

Question

How to take the variance of a second order expansion? $\text{Var}\left[aX+bY+cXY+mX^2+nY^2\right]$

Let say we have 5 real-valued constant parameters $\{a,\ b,\ c,\ m,\ n\}$, and two random variables $X$ and $Y$ such their means values are $\mu_x$ and $\mu_y$, their variances are $\sigma_x^2$ and $\sigma_y^2$, and they are correlated with coefficient $\rho_{xy}$.

I want to know how to calculate the following variance and which value it has based on the previous parameters: $$\text{Var}\left[aX+bY+cXY+mX^2+nY^2\right]$$

1. Could this be done if "nothing" is known about the probability distribution of each random variable (maybe different distributions for each one)? (but considering that it is true that each one has finite parameters $\mu_i$, $\sigma_i^2$, and $\rho_{ij}$).
2. If point (1) is not possible: Could this be answer if $X$ and $Y$ are normally distributed? (I found the Isserlis' theorem and Stein's lemma, but I got stuck anyway).

I aiming to find a formula I could directly use, so please explain if your results are based in any assumptions (unfortunately, independence is not applicable for the result I am looking for).

Answer 1

Point (1) is impossible, because the given variance can assume any value (or may fail to exist) if there is no further constraints on the joint distribution of $X$ and $Y$.

As for (2), the answer is affirmative if $X$ and $Y$ are jointly normal (i.e., multivariate normal). In this case, it is easy to check that

$$ (X, Y) \quad \stackrel{d}= \quad \bigl(\mu_X + \sigma_X Z, \mu_Y + \sigma_Y(\rho_{XY} Z + \sqrt{1-\smash[b]{\rho_{XY}^2}} W) \bigr)$$

for i.i.d. standard normal $Z$ and $W$. Using this, we can in principle compute the variance of OP's variable

$$S = aX + bY + cXY + mX^2 + nY^2.$$

However, the answer turns out to be nasty. For simplicity, let

\begin{align*} \mu_X^{(2)} &= \mathbf{E}[X^2] = \mu_X^2 + \sigma_X^2, \\ \mu_Y^{(2)} &= \mathbf{E}[Y^2] = \mu_Y^2 + \sigma_Y^2, \\ \sigma_{XY}^{(2)} &= \mathbf{Cov}(X, Y) = \rho_{XY}\sigma_X\sigma_Y. \end{align*}

Then $\mathbf{Var}(S)$ can be expanded as:

(Forgive my laziness for choosing to upload a formula instead of typing it out. This is too long, though...)

Answer 2

For the first point, not, as it was explained in the comments. For the 2nd point, if we can assume normality, it is tedious and perhaps not so easy: $$ \text{Var}\left[aX+bY+cXY+mX^2+nY^2\right] = \mathbb{E}\left( aX+bY+cXY+mX^2+nY^2 - \mu \right)^2, $$ where $\mu = \mathbb{E}\left( aX+bY+cXY+mX^2+nY^2\right) = a\mathbb{E}(X) + b\mathbb{E}(Y) + c\mathbb{E}(XY) + m\mathbb{E}(X^2) + n\mathbb{E}(Y^2)$.

So now expand the parentheses to find the formula: $$ \text{Var}\left[aX+bY+cXY+mX^2+nY^2\right] = a^2\mathbb{E}(X^2) + b^2\mathbb{E}(Y^2) + c^2\mathbb{E}(X^2Y^2) + \dots - \mu^2. $$ I'm not going to bother with this simple expansion. I'm sure you can derive that yourself. You will then need to know the expectations of products of normal distributions. The product of two normal r.v.s is known and can be expressed as a chi-squared distributed r.v.. Or find it explicitly here. For the product of three or four normal r.v., it's much harder to get something explicit. You can try to read this paper.

Addendum: I'm not assuming joint normality like the other answer.