Is the prox-residual monotone? $\newcommand{\scp}[2]{\langle #1,#2\rangle}\newcommand{\id}{\mathrm{Id}}$
Let $f$ and $g$ be two proper, convex and lower semi-continuous functions (on a Hilbert space $X$ or $X=\mathbb{R}^n$) and let $g$ be continuously differentiable. Consequently, the subdifferential $\partial f$ and the gradient $\nabla g$ are monotone operators, i.e.
$$
\scp{\nabla g(x)-\nabla g(y)}{x-y}\geq 0
$$
for all $x,y$ and for $u\in\partial f(x)$ and $v\in\partial f(y)$
$$
\scp{u-v}{x-y}\geq 0.
$$
My question:

Is the operator $x\mapsto x - (\id + \gamma\partial f)^{-1}(x-\gamma\nabla g(x))$ monotone for $\gamma\geq 0$?

Note that I ask for all values $\gamma\geq 0$ specifically for large values.
An equivalent question is

Does for $\gamma\geq 0$ it holds that
  $$
\scp{(\id + \gamma\partial f)^{-1}(x-\gamma\nabla g(x)) - (\id + \gamma\partial f)^{-1}(y-\gamma\nabla g(y))}{x-y}\leq \|x-y\|^2
$$

Thoughts:
For $f=0$ and $g=0$ it's clear (for $g=0$ this follows since the "proximal operator" $x\mapsto (\id + \gamma\partial f)^{-1}(x)$ is non-expansive, for $f=0$, the monotonicity of $\gamma\nabla g$ works in the right direction and does the trick). 
Also for small $\gamma$ (smaller that $2/L$ if $L$ is the Lipschitz constant of $\nabla g$, if is has one) the thing is clear as then the mapping $x\mapsto (\id + \gamma\partial f)^{-1}(x - \gamma\nabla g(x))$ is again non-expasive (and used as iteration in the so-called proximal gradient method). However, for large $\gamma$ I could not make use of any of these observation since using Cauchy-Schwarz for the inner product ruins the estimate then.
Moreover all examples I tried numerically (in various dimensions and for various functions $f$ and $g$) suggested that the claim holds.
Intuitively, all these together makes me think that the answer to my questions is yes, but I failed to prove it. Also all inequalities in Bauschke/Combettes "Convex analysis and Monotone Operator Theory in Hilbert spaces" I found were not helpful. 
 A: This inequality is not true. Here is a counterexample:
Let ${\bf x},{\bf y} \in \mathbb{R}^2$ with ${\bf y} = 0$ and ${\bf x} = (1,\epsilon)$ for $\epsilon$ small positive. Let $g(x,y) = Cy^2$ with $C(\epsilon)$ very large, chosen say so that
$${\bf x} - \nabla g({\bf x}) = (1, -M)$$
for some $M$ large.
We now construct $f$ convex so that $\nabla(|x|^2/2 + f)(2,-10) = (1,-M)$ and $\nabla f(0) = 0$, which would violate the desired inequality. We rewrite this as $\nabla f(2,-10) = (-1,10-M)$, $\nabla f(0) = 0$. Take
$$f(x,y) = \frac{1}{2}(C_1y^2 + C_2(x+y)^2)$$
for $C_1,C_2$ to be chosen momentarily. Then
$$\nabla f(2,-10) = (-8C_2, -10C_1 - 8C_2), \quad \nabla f(0) = 0.$$
Taking $C_2 = 1/8$ and $C_1$ so that the second component is $10-M$ (if $M$ is large then $C_1$ can be chosen positive, making $f$ convex), we are done.
The constants in this example are somewhat arbitrary; geometrically, the point is that we can make $g$ "very monotone" in one direction so that the non-expansivity of $(Id + \nabla f)^{-1}$ can't save us.
A: Although not monotone at the operator level (as suggested by C. Mooney's proof), the monotonicity of prox-residual norms is known (probably you are already aware of it). 
Let $P_\eta^g$ denote the prox-operator for $g$ with index $\eta$, i.e.,
\begin{equation*}
 P_\eta^g(y) := \text{argmin}_{x} \tfrac{1}{2\eta}\|x-y\|^2 + g(x).
\end{equation*}
Let $y, z \in \mathbb{R}^n$ and $\eta > 0$. Define the functions
\begin{eqnarray*}
p_g(\eta) &:=& \tfrac1\eta\|P_\eta^g(y-\eta z) - y\|\\
q_g(\eta) &:=& \|P_\eta^g(y-\eta z) - y\|.
\end{eqnarray*}

Theorem. $p_g$ is a decreasing function of $\eta$ and $q_g$ is an increasing function of $\eta$.

