Can i include the product of two random variables? Or do I risk collinearity? I have a model in which I want to predict Y.
My regressors X, are x1 and x2.
For some reason I believe that it would also be useful to include into the model:


*

*x3 = x1 * x2

*x4 = x1 / x2
Can I use a regressors x1, x2, x3 and x4 altogether or do I risk perfect collinearity problem.
I know for instance that using x5 = x1 + x2 would yield perfect collinearity and hence a completely useless regressor.
 A: You won't have perfect collinearity (as per your question), but you do risk multicollinearity issues with your two additional regressors.
While they're not algebraically linear combinations of the two predictors, it can be the case that these variables (x1-x4) in a particular sample might lay close to to a linear subspace - with the typical consequences of near-multicollinearity.
For example, if the two original variates both have very small coefficients of variation then their product can be quite closely related to their sum (or some other linear combination if they're dissimilar in size). This can happen even if the original variables are not highly correlated.
Similarly, high correlation can happen with the ratio and the difference.
A: Operating under your stated assumption that $x_3=x_1x_2$ and $x_4=x_1/x_2$ need to be entertained as possible explanatory variables in a model of a response $Y$ (and therefore not summarily dropped because they might be a little inconvenient), it can be helpful to consider alternative ways of expressing this model.
As stated, the model is of the form
$$Y \sim F(x_1, x_2, x_3, x_4; \theta) = F(x_1, x_2, x_1x_2, x_1/x_2; \theta)$$
for a given distribution family $F$ involving unknown parameters $\theta$ to be determined.  For instance, a linear regression model would involve a five-dimensional parameter $\theta = (\beta_0, \beta_1, \beta_2, \beta_3, \beta_4)$ in the form
$$E[Y] = \beta_0 + \beta_1 x_1 + \beta_2 x_2 + \beta_3 x_3 + \beta_4 x_4 = \beta_0 + \beta_1 x_1 + \beta_2 x_2 + \beta_3 x_1x_2 + \beta_4 x_1/x_2.$$
For simplicity of exposition, let's analyse the linear regression model: it will be clear how the analysis extends to other models.
One way is to restate the model in terms of $x_4$ and $x_2,$ which algebraically imply $x_1=x_2x_4$ and $x_3=x_2^2x_4:$
$$E[Y] = \beta_0 + \beta_1 x_4x_2 + \beta_2 x_2 + \beta_3 x_4x_2^2 + \beta_4 x_4 =  \beta_0 + \beta_2 x_2 + \beta_4 x_4 + x_4\left(\beta_1 x_2 + \beta_3 x_2^2\right).$$
The last term would ordinarily be characterized as an interaction between $x_4$ and a quadratic function of $x_2.$  Since, except in very special circumstances, interactions should be included only when their component terms are included, this suggests you ought to extend to model to include an $x_2^2$ term.  It would have the form
$$E[Y] =  \beta_0 + \left(\beta_2 x_2 + \beta_5 x_2^2\right) + \beta_4 x_4 + x_4\left(\beta_1 x_2 + \beta_3 x_2^2\right).$$
That is a model involving (a) $x_4$ and (b) the simplest possible quadratic spline of $x_2.$ Such models are common: the quadratic terms allow for some amount of nonlinear response in $x_2$ and the interaction allows for the response to change with different values of $x_4$ in a controlled way.
These simple algebraic manipulations demonstrate that the proposed model is not at all unusual.  They reframe it in terms of standard, well-understood concepts.

There remains the question of collinearity.  That collinearity could be a problem is demonstrated by the case where both $x_1$ and $x_2$ are binary variables coded as $\pm 1.$  In this case, $x_1/x_2$ and $x_1x_2$ are always equal (not just collinear).  
On the other hand, that collinearity might not be much of a problem can be demonstrated by exhibiting some sample data with relatively little collinearity.  We would want $x_2$ to be orthogonal to $x_2^2,$ of course, and then everything will be ok provided the interactions don't introduce collinearity.  Unfortunately, $x_4$ and $x_4x_2^2$ are likely to be positively correlated.  But by how much?
Consider the data $x_2 = (-1,0,1,\, -1,0,1,\, -1,0,1)$ and $x_4 = (-1,\sqrt{3},-1,\,0,0,0,\,1,-\sqrt{3},1).$  The covariance matrix of the columns $(x_2, x_2^2, x_4x_2, x_4, x_4x_2^2)$ is 
$$\pmatrix{3&0&0&0&0 \\ 0 & 1 & 0&0&0 \\ 0&0&2&0&0 \\ 0&0&0&5&2 \\ 0&0&0&2&2}/4.$$
It is nearly orthogonal, with correlation only between the last two variables (as expected).  (Notice that introducing $x_2^2$ has not changed anything, because this variable is orthogonal to all the others.) The ratio of the largest to the smallest eigenvalue (its condition number) is $6.$  This is not beautiful, but it's not bad, either.  One could easily obtain reliable coefficient estimates with such explanatory variables.
If you don't have the luxury of choosing the values of $x_2$ and $x_4$ to arrange such near-orthogonality, then you will simply have to proceed as anyone would always do in such cases: investigate the data you have and deal with any collinearity in the usual ways (which would include ignoring it; dropping variables based on scientific considerations; selecting some principal components; using a Lasso; and so on).
A: No, you don’t risk collinearity because $x_i$ are not linearly dependent in general, i.e. the below equation has just one solution holding for all possible $x_i$:
$$a_1x_1+a_2x_2+a_3x_3+a_4x_4=0$$
And that is $a_i=0$.
In $x_5=x_1+x_2$ case, the following equation has non-zero solutions such that $a_1=a_2=-a_5$: $$a_1x_1+a_2x_2+a_5x_5=0$$
