Gravitational potential energy of any spherical distribution
Start by writing that: $$\tag1 M=\int^R_0 \rho(r)4\pi r^2 dr = \int^R_0 \frac{dM}{dr}dr.$$ From (1), we can identify that: $$\tag2 \frac{dM}{dr}dr = 4\pi r^2 \rho(r)dr.$$ Now let's return to your equation (1): $$U = - \int_0^R \frac{GM(r)}{r} \, \rho(r) \, 4 \pi r^2 dr,$$ using equation (2), this becomes: $$U=-G\int_0^R \frac{M(r)}{r} \frac{dM}{dr}dr.$$
Let's concentrate on just the integral part, and call it $I$ for concreteness, $$I=\int_0^R \tag{3}\frac{M(r)}{r} \frac{dM}{dr}dr.$$ Aha! This looks like something we can integrate by parts. Let's set $u=\frac{M(r)}{r}$, and $v'=\frac{dM}{dr}$. Then $$u'=-\frac{M(r)}{r^2} + \frac{1}{r}\frac{dM}{dr}.$$ Continuing with IBP, we get: $$I=\Bigg[\frac{M(r)}{r}M(r)\Bigg]^R_0-\int^R_0 \Big(-\frac{M(r)}{r^2} + \frac{1}{r}\frac{dM}{dr}\Big)M(r) dr.$$ Let's look at solely the integral part, calling it $J$: $$J=-\int^R_0 \frac{M^2(r)}{r^2} dr + \underbrace{\int^R_0 \frac{M(r)}{r}\frac{dM}{dr} dr}_{\textrm{we know this integral, it's }I}.$$ So now we have that: $$I=\Bigg[\frac{M^2(r)}{r}\Bigg]^R_0 + \int^R_0 \frac{M^2(r)}{r^2} dr - I,$$ and so: $$2I= \Bigg[\frac{M^2(r)}{r}\Bigg]^R_0 + \int^R_0 \Bigg(\frac{M(r)}{r}\Bigg)^2 dr.$$ For the case of non-pathological $M$, that is $M(0)=0$, we get: $$I=\frac{M^2}{2R} + \frac{1}{2}\int^R_0 \Bigg(\frac{M(r)}{r}\Bigg)^2 dr.$$
Hence we can identify from this that $U$ has the limit of $-\frac{GM^2}{2R}$, and from that we subtract the integral of something positive definite.
I believe this proves (to physicist standards) your conjecture about $k$.
Request of simplifying the process:
Let's just consider with no dimensions the integral: $$K=\int^X_0 \frac{f(x) f'(x)}{x} dx.$$
Take: $$g(x)=\frac{f^2(x)}{2x},$$ then taking the derivative: $$g'(x) = \frac{f(x) f'(x)}{x}-\frac{f^2(x)}{2 x^2}.$$ We can recognise our original integral in this mix, and so integrating $g'(x)$ will give us our integral, less $\int^X_0 \frac{f^2(x)}{2 x^2} dx$, and so we can write: $$K=\big[g(x)\big]^X_0 + \frac{1}{2}\int^X_0 \frac{f^2(x)}{x^2}dx= \Bigg[\frac{f^2(x)}{2 x}\Bigg]^X_0 + \frac{1}{2}\int^X_0 \bigg( \frac{f(x)}{x} \bigg)^2 dx.$$ Giving us the result we wanted. As IBP is always the inverse of the product rule, with a little bit of thought you can generally figure out the original function. However if I'd just posted this last 4-step proof, I think you'd have just thought I was a wizard...
Your claim, $k\in\bigg[\dfrac{1}{2},\infty\bigg)$, is true.
I will provide a hand-wavy proof that physicists might accept as a proof but mathematicians will surely ridicule as a conjecture. So, here it goes:
If invert the sign of your integral (and we will set $4\pi G=1$) then it reads: $$\displaystyle\int\dfrac{M_{\text{in}}(r)}{r}\rho(r) r^2 dr \tag{1}$$ This is identical to the expression (provided that we have set $4\pi\epsilon_0=1)$ for the self-potential energy of a spherical charge distribution if we identify $M_{\text{in}}(r)$ with the charge within radius $r$ and $\rho(r)$ with the charge density at radius $r$. Now, we know that if we put a certain charge on a metallic spherical ball then, in its equilibrium state, the charge distribution will be so that it minimizes the self-potential energy that a charge distribution of total charge $M$ can have when distributed within a radius $R$. But we know, according to Gauss's law, that the charge gets distributed in the form of a spherical shell on the surface of the spherical ball. Thus, the minimum value of $(1)$ will be when all the charge is at $R=1$. We know, as you calculated, the self potential energy in this case becomes $\dfrac{1}2\dfrac{M^2}R$. Of course, the maximum value $(1)$ can have is indefinitely large corresponding to larger and larger fraction of charge being put nearer and nearer to the center.
In the gravitational case, the additional minus sign makes you read minimum of the above explanation as maximum and vice-versa. But since you have kept the negative sign outside of the definition of $k$ (i.e. you multiply $k$ with a minus to get the final $U$), you can very directly see that $k\in\bigg[\dfrac{1}{2},\infty\bigg)$ holds.