I have the following problem. Let $\mu$ and $\nu$ be two probability measures on $\mathbb{R}$ and let $f:\mathbb{R}\times \mathbb{R}\to \mathbb{R}_+$ be a non-negative continuous function. Assume that we have two non-negative functions $g,h:\mathbb{R}\to \mathbb{R}_+\cup \{+\infty\}$ such that $g$ is in $L^1(\mathbb{R},\mu)$ and $h$ is in $L^1(\mathbb{R},\nu)$ verifying $$f(x,y) \leq g(x) + h(y) $$ for all $x,y$ in $\mathbb{R}$. Note that the sum $g(x) + h(y)$ is always well defined.
Even if I'm not completely convinced, it seems intuitive to me that, since $f$ is finite everywhere, it should be possible to take $g$ and $h$ to be finite everywhere. Nevertheless, I could not find a proof nor a counterexample. Any suggestions?
This is too long for a comment and only a hint/suggestion, in other words not a solution.
We assume in what follows that $\mu$ and $\nu$ are not fixed or given in advance, and that only $f: \mathbb{R} \to \mathbb{R}^+$ is fixed (I know that the domain is wrong -- hopefully this generalizes to that case). Then we show how, under nice circumstances, we can find probability measures $\mu$ and $\nu$ and functions $g$ and $h$ which may satisfy the required conditions.
Define $f$ to be the cumulative distribution function of some signed measure $\omega$ on $\mathbb{R}$, for all $x \in \mathbb{R}$, $$\omega(-\infty,x]=f(x) $$ Since the sets of the form $(-\infty,x]$ generate the Borel $\sigma-$algebra on $\mathbb{R}$, this is enough to specify a signed measure $\omega$ on the measurable space $\mathbb{R}$ with the Borel $\sigma$-algebra. The Lebesgue $\sigma-$algebra is just the completion of the Borel $\sigma-$algebra hence presents no new complications.
Note that $\omega$ is not necessarily absolutely continuous with respect to Lebesgue measure; that follows only if $f$ is an absolutely continuous function. For example, $f$ could be the cumulative distribution function of the Cantor distribution/measure. This still wouldn't be a problem though from the perspective of the Lebesgue decomposition theorem.
Writing $f$ as the cumulative distribution of some signed measure $\omega$ is possible because of the Lebesgue decomposition theorem and the fact that $f$ is a continuous function.
In general we have to assume that $\omega$ is a signed measure, as opposed to a non-negative measure, because we have no assumption that $f$ is monotone increasing. Even being non-negative, $f$ can still decrease, just not indefinitely.
The assumption that $f$ is finite everywhere means that $\omega$ is $\sigma-$finite, even if $f$ is unbounded as $x \to \infty$ or $x \to -\infty$. Therefore the Hahn decomposition theorem applies.
Then by the Hahn decomposition theorem, there exist two measures $\mu$ and $\nu$ on $\mathbb{R}$ such that $$\omega = \mu^* - \nu^*$$
If $\mu^*$ is a finite measure, then its cumulative distribution function is well-defined. Likewise, since $f$ is non-negative, and thus cannot decrease indefinitely, if $\mu^*$ is a finite measure, than $\nu^*$ must also be a finite measure, because we must have for all $x$ that $\nu^*(-\infty,x] \le \mu^*(-\infty,x]$.
(Since $\mu^*(-\infty,x] - \nu^*(-\infty, x] = \omega(-\infty, x] =f(x) \ge 0$.)
However, if $\mu^*$ is not a finite measure, then it cannot be normalized to a probability measure. Moreover, $\mu^*$ does not necessarily have to (although it can) have a well-defined cumulative distribution function if it is not a finite measure.
Likewise, if $\mu^*$ is not a finite measure, then $\nu^*$ can also be (but does not have to be) an infinite measure, provided that $\nu^*[a,b] \le \mu^*[a,b]$ for all $a<b$. Then $\nu^*$ also could not be normalized to a probability measure, and it would not need to have a well-defined cumulative distribution function.
That being said, if $\mu^*$ is an infinite measure, it doesn't matter what $\nu^*$ is, because we are already stuck. And again, if $\mu^*$ is a finite measure, then $\nu^*$ must be also.
So in order to proceed, we must assume that $\mu^*(\mathbb{R}) < \infty \implies \nu^*(\mathbb{R}) < \infty$. So now define probability distributions on $\mathbb{R}$ by: $$\mu = \frac{1}{\mu^*(\mathbb{R})} \mu^*$$ $$\nu = \frac{1}{\nu^*(\mathbb{R})} \nu^*$$ If we define $\tilde{g}$ to be the cumulative distribution function of $\mu$ and $\tilde{h}$ to be the cumulative distribution function of $\nu$, then we have for $$g:= \mu^*(\mathbb{R}) \tilde{g}$$ $$h:=-\nu^*(\mathbb{R}) \tilde{h}$$ that $$f=g+h$$ and presumably that $g \in L^1(\mathbb{R}, \mu)$ and $h \in L^1(\mathbb{R},\nu)$ (although I am not sure how to evaluate the corresponding integrals, so I can neither confirm nor deny that they are indeed finite).