Approximation for Newton's formula for gravity as modified by the expansion of the universe? Newton's formula for the gravitational force between two objects is
$$F = -\frac {Gm_1m_2}{r^2}.$$
Assuming the force is due to the exchange of gravitons and that gravitons are affected by the expansion of the universe in the same way as photons, can I simply modify the formula so it behaves in the same way as radiation pressure?
$$F = -\left(1-\frac {Hr} c\right)\frac {Gm_1m_2}{r^2}$$
where $H$ is the Hubble constant and $c$ is the speed of light
I suppose it also assumes that the force is small and that the objects are not gravitationally bound.  Is this a reasonable first order approximation?
[Edit 3] I appreciate the answers, but I don't feel comfortable with them. A uniform constant expansion of the universe gives rise to an acceleration of $H^2r$ between any two points separated by $r$. This is because $r$ increases with time and so the recession velocity Hr also increases with time. But this does not cause any acceleration with respect to the CMB, for example, and is not a function of either mass and cannot be included as a correction to Newtons formula for force.
If we define force as mass times acceleration with respect to the CMB, then I'd like to suggest
$$F = - \left (1-\frac {Hr} c  \right )^4 \frac {Gm_1m_2}{r^2} $$
with the first order correction being approximately $$ 4 \frac H c \frac {Gm_1m_2} r $$
This expression is simply by analogy with the radiation pressure from a distant black body.
 A: starting with this metric
\begin{align*}
&\mathbf{G}=\left[ \begin {array}{cccc} B \left( r \right) {c}^{2}&0&0&0
\\  0&-A \left( r \right) &0&0\\  0&0
&-{r}^{2}&0\\  0&0&0&-{r}^{2} \left( \sin \left(
\theta \right)  \right) ^{2}\end {array} \right] 
\end{align*}
you obtain from the geodetics equations and the first integrals
\begin{align*}
 &r^2\,\frac{d\phi}{d\tau}=l=\text{constant}\\
 &B(r)\frac{dt}{d\tau}=\frac{1}{c}\,f=\text{constant}
\end{align*}
that (conservation of the energy)
\begin{align*}
  &A(r)\left(\frac{dr}{d\tau}\right)^2+\frac{l^2}{r^2}-\frac{f^2}{B(r)}=-\epsilon=\text{constant}\tag 1
\end{align*}
with $~\epsilon=c^2~$ (from  $~ds^2~$ ) and with  the Schwarzschild metric
\begin{align*}
 &B(r)=\frac{1}{A(r)}= 1-\frac{2\,a}{r}\quad,a=\frac{G\,M}{c^2}\quad\Rightarrow
\end{align*}
equation (1)
\begin{align*}
 &\frac{1}{2}\left(\frac{dr}{d\tau}\right)^2+U(r)=\frac{f^2-\epsilon}{2}=\text{constant}\\
 &U(r)=-\frac{G\,M}{r}+\frac{l^2}{2\,r^2}-\frac{G\,M\,l^2}{c^2\,r^3}
\end{align*}
thus the Newton correction  is $~-\frac{G\,M\,l^2}{c^2\,r^3}$

with the cosmological constant $~ \Lambda$
$$B(r)\mapsto B(r)-\Lambda\,\frac{r^2}{3}$$
thus the Newton correction is now
$$-\frac{G\,M\,l^2}{c^2\,r^3}-\Lambda\,\frac{r^2}{6}$$
Force
with $~l=\frac Lm~$
$$F=-m\,\frac{dU}{d\tau}=
{\frac {{L}^{2}}{m{r}^{3}}}-3\,{\frac {G\,M{L}^{2}}{m{c}^{2}{r}^{4}}}-{
\frac {m\,M\,G}{{r}^{2}}}+\frac 13\,m\,r\Lambda
$$
De-Sitter  solution
with $~H^2=\frac {\Lambda}{ 3}\quad\Rightarrow$
$$F=
{\frac {{L}^{2}}{m{r}^{3}}}-3\,{\frac {G\,M{L}^{2}}{m{c}^{2}{r}^{4}}}\underbrace{-{
\frac {m\,M\,G}{{r}^{2}}}+m\,r\,H^2}_{=-(1-x)\frac{m\,M\,G}{r^2}}
$$
thus
$$x=\frac{r^3\,H^2}{M\,G}$$
Friedmann equation for materiel dominate cosmos  ($~c=1~$)
\begin{align*}
& \dot{r}^2-\frac{K_m}{r}-\frac{1}{3}\Lambda\,r^2=-k\tag 1
\end{align*}
with
\begin{align*}
 &K_m=\frac{8\,\pi\,G}{3}\rho_m\,r^3=\text{constant}
\end{align*}
multiply equation (1) with $~\frac M2~$ you obtain the energy
\begin{align*}
 &\frac{M}{2}\dot{r}^2\underbrace{-\frac{G\,M^2}{r}-\frac{M}{2}\,\frac{1}{3}\Lambda\,r^2}_{U(r)}=-\frac{M}{2}\,k=\text{constant}
\end{align*}
where
\begin{align*}
 &M=\frac{4\,\pi}{3}\rho_m\,r^3
\end{align*}
thus the force is (with $~\Lambda=3\,H^2~$)
\begin{align*}
  &F=-\frac{G\,M^2}{r^2}+M\,H^2\,r
\end{align*}
A: The gravitational force doesn't redshift. It Lorentz contracts:

This image shows the electric field, but the (relativistic) gravitational field behaves similarly. The gravitational field of an object moving away from you is effectively weaker, but so is that of an object moving toward you,
Perhaps you could try to capture that in a version of the force law that includes a dependence on velocity as well as position, but if you want to go that far in the pursuit of correctness, then I think you should also use the retarded position and velocity. If you do everything correctly, you will end up with linearized GR.
$H$ shouldn't be in the formula, because the strength of gravity doesn't depend on it. $Hr$ may be a good approximation to an object's recessional velocity, but since you're writing a per-object force law (where $r$ is the distance between two particular objects), you should use the actual velocity of those objects, not an approximation to it.
I think the best quasi-Newtonian force law is the one in the Q&A that J.G. linked in a comment. You can derive the second Friedmann equation with $p=0$ from it, and since $p\approx 0$ for most of the universe's history, the Newtonian cosmology is a pretty good approximation to ΛCDM.
A: There is an exact solution where (in $G=c=1$ units, which everyone should be using for doing general relativity) you replace $1-2M/r$ with $1-2M/r - \frac{1}{3}\Lambda r^{2}$, which will give you $T_{ab} = \Lambda g_{ab}$ (you can verify by computing the Ricci curvature of the above metric directly), and so we can think of as "asymptotically (anti) de Sitter Schwarzschild".
Since for de Sitter in 3+1 dimensions, $a(t)= a_{0}e^{t\sqrt{\Lambda/3} }$, we have $H = \sqrt{\Lambda/3}$, and at least for this case, a better ansatz would be to generalize Newton by $F_{{\rm on\;} m_{1}} = Gm_{1}\frac{m_{2}}{r^{2}} - 2m_{1}H^{2} r$ rather than the multiplicative rule you took.
