Interaction term in multivariate polynomial regression I'm looking for answer for the question about multivariate polynomial regression.
I can't find a clear explanation of when an interaction term is necessary.
Some sources say that the estimated model of a complete second degree polynomial regression model in two variables $x_{1}$, $x_{2}$ may be expressed as
$$\hat{y} = b_{0} + b_{1}x_{1} + b_{2}x_{2} + b_{3}x_{1}^{2} + b_{4}x_{2}^{2} + + b_{5}x_{1}x_{2}$$
and others don't consider interaction term $x_{1}x_{2}$...
And how does it look when I have for example $d=4$ and $p=3$, a third degree polynomial of four independent variables?
 A: One could use the model with or without interaction terms, depending on whether one expected them to be useful. In practice, I would expect any description of a polynomial regression model to be clear about whether or not interaction terms were included. If no mention is made of them, I would assume not.
A: One consideration in such problems is the interpretability of effects. The interpretation of higher order terms as well as their interactions is somewhat involved, but important to understand. These are discussed in other CV posts, so do check them. 
The consensus is that, with higher order terms as well as their interactions, you may fit interactions as needed provided that the order does not exceed the order of the main effects. This preserves the interpretation of the effect of product terms as a difference in differences. So in your example, it makes sense to include the $x_1 x_2$ product term since you have their respective main effects $x_1$ and $x_2$. Furthermore, you could also include $x_1 x_2^2$, $x_1^2 x_2$, and/or $x_1^2 x_2^2$ since they are all nested within the lower level terms. One would not, however, adjust for (say) $x_1^3 x_2$  without first adjusting for the $x_1^3$ main effect.
The number of possible product terms in a model controlling for polynomial orders quickly becomes untenable. A "semi-saturated" model for the outcome is quickly approached with relatively few effects. For instance, if $d$ is the polynomial order of $x_1$ and $q$ is the polynomial order of $x_2$ then there are $d + q$ main effects and $d \times q$ possible product terms. Fitting semi-saturated models has a purpose in some estimation problems (like robust variance estimation), but you quickly approach an issue of overfitting in most settings. Using an information criterion like the BIC can help to select a model with a lower number of overall effects, but high predictive accuracy and external validity.
