It's not hard to see that $H$ must be of the form $\alpha ux + \beta uz + \gamma yx + \delta yz$ by $\mathbb{R}$-linearity of $F$ (see Jan-Cristoph Schlage-Puchta's answer for a fuller explanation).

By using symmetry of multiplication, we can assume without loss of generality that $\beta = \gamma$.

We have the equation $$F(P \times Q) = H(F(P), F(Q), P, Q) = \alpha F(P) F(Q) + \beta F(P) Q + \gamma P F(Q) + \delta PQ$$

Setting $P = Q = 1$, we get $F(1) = \alpha F(1)^2 + \beta F(1) + \gamma F(1) + \delta$. Degree-checking tells us that if $F(1)$ isn't constant, then $\alpha = 0$. More careful checking tells us that $\beta + \gamma = 1, \delta = 0$. This leads to the solution $F_{M, R}(P) = P \times R$, where $R = F(1)$. That answer aside, we assume $F(1)$ is constant.

I'm going to take a detour for a moment, to talk about an action of the two-dimensional nonabelian Lie group on the space of pairs $(F, H)$ that satisfy $F(P \times Q) = H(F(P), F(Q), P, Q)$. This action comes from the map $g_{(a, b)} F = aF + b \operatorname{Id}$. We want the equation to remain the same, which leads us to:

$$(gH)(u, x, y, z) = aH(\frac{u - by}{a}, \frac{x - bz}{a}, y, z) + b yz$$

$$ = \frac{\alpha}{a} ux + (\beta - \frac{b}{a} \alpha) uz + (\gamma - \frac{b}{a} \alpha) yx + (a \delta - b \beta - b \gamma + \frac{b^2}{a} \alpha + b) yz$$

Back from the detour: we now set just $Q = 1$, getting $F(P) = \alpha F(P) F(1) + \beta F(P) + \gamma P F(1) + \delta P$. Rearranging, we get $(1 - \alpha F(1) - \beta) F(P) = (\gamma F(1) + \delta) P$. The obvious solutions have $F(P) = \lambda P$ for some $\lambda \in \mathbb{R}$, which have been discussed in other answers; otherwise, we must have $1 - \alpha F(1) - \beta = \gamma F(1) + \delta = 0$.

We now split into two cases: either $\alpha = 0$, or $\alpha \neq 0$. If $\alpha = 0$, then $\beta = \gamma = 1$. Using the action of the group, we can reduce to the case that $\delta = 0$, implying that $F(P \times Q) = F(P) Q + P F(Q)$, i.e. that $F$ is a derivation. This can only happen when $F(P) = R \times \partial P$ for some polynomial $R$. Undoing the group action, we get $F_{D, R, \lambda}(P) = \lambda P + R \times \partial P$.

If $\alpha \neq 0$, then we can use the action of the group to reduce to the case $\alpha = 1, \beta = \gamma = 0$. But then $\delta = 0$, giving us the equation $F(P \times Q) = F(P) F(Q)$ - in other words, $F$ is a homomorphism. Homomorphisms from $\mathbb{R}[X]$ are compositions: $F(P) = P \circ R$ for some polynomial $R$. Undoing the group action, we get $F_{H, R, \lambda} = \lambda P + P \circ R$.

So all solutions $(F, H)$ are of the forms:

1) $F(P) = \lambda P, \frac{1}{\lambda} \alpha + \beta + \gamma + \lambda \delta = 1$

2) $F(P) = \lambda P + R \times \partial P, \alpha = 0, \beta = \gamma = 1, \delta = \lambda$

3) $F(P) = \lambda P + c (P \circ R)$ with coefficients that aren't difficult to determine using the group action.

4) $F(P) = P \times R, \alpha = 0, \beta = \gamma = \frac{1}{2}, \delta = 0$