To fix the value of the coupling coefficient $κ$ that appears in Einstein's equation $G_{μν} = κ T_{μν}$ and connect it to the constant $G$ that appears in Newton's law of gravity, you need the correspondence limit with Newtonian gravity. In that way, you ensure backwards compatibility.
So, in that sense, you need Newton's law of gravity ... but just to pin down the coefficient.
It is a value that most on-line cites state wrong, because of the propagation of widespread errors, and a value which the literature itself gets confused for the same reason. I've seen 3 different values published in the literature. There will be more to say on that below.
The simplest way to establish a correspondence is to solve the vacuum equation $R_{μν} = 0$ for a spherically symmetric source - to represent the exterior of a spherical clump of matter, such as a planet, and then find what its Newtonian analogue would be. There is no coupling coefficient in the vacuum equation, so you'll be able to establish the geometry first, before drawing correspondences, without worrying about what $κ$ should be.
In a way, that was the idea you already conveyed.
The solution you find is the Schwarzschild solution and it has a form given by the line element in spherical coordinates
$$\frac{dr^2}{1 - k/r} + r^2\left(dθ^2 + (\sin θ dφ)^2\right) - c^2 dt^2 \left(1 - \frac{k}{r}\right). \label{Rel0}\tag{Rel0}$$
On the exterior of the body, $k$ is constant. Inside the body, $k$ varies with $r$, from $k = 0$ at $r = 0$, up to its exterior value at $r = R > 0$, where $R$ is the radius of the body. It is the exterior part we are interested in, here, and the exterior value for $k$.
Now, we compare it with the Newtonian version and draw correspondences, establish what the model parameter $k$ ought to be, and then proceed from there to then find what the coefficient $κ$ should be. But, there might be a minor hang-up at this point: what's the Newtonian version of this?! Well, there is an answer - a simple one - but not quite what anyone might be expecting.
You see, Einstein's equation is really a transport law for fluid dynamics with the gravitational effect of curved space-time geometry built into it. The left hand side of the equation $G_{μν}$ - Einstein's tensor - is geometrically-based. On the right, you have $T_{μν}$, which contains the properties of matter ... in continuum form: as densities and fluxes. That includes mass density, momentum density, energy density and their respective fluxes, as well as stress and pressure in the case of momentum density.
Actually, it doesn't ...
... and that's one of the errors that has been widely propagated in the literature and is the root of a lot of confusion. The actual quantity that contains the various densities and fluxes (as the term "density" already clearly indicates) is not a tensor at all, but a tensor density. Here, it's the stress tensor density $𝔗^ρ_ν$, which the stress tensor is related to by: $\sqrt{|g|} T_{μν} = g_{μρ} 𝔗^ρ_ν$, where $g_{μν}$ are the components of the spacetime metric and $g$ is the determinant of the matrix formed by them.
But that's a digression from the main point here (and a digression we'll come back to in a moment). The main issue is this: in Relativity, mass and kinetic energy are combined into "total energy". Call the respective quantities $m$, $H$ and $E$. What one finds is that there is a contribution to $m$ that comes from the kinetic energy $H$ that yields a "moving mass" $M = m + H/c^2$. The total energy for a body is then given by That Equation $E = Mc^2$; i.e. $E = H + mc^2$.
So, what was 5 now becomes 4: the total energy $E$ and the three components of the momentum $𝐩 = \left(p_x, p_y, p_z\right)$, instead of $M$, $H$ and $𝐩$. They are combined into a single vector, with $p_0 = -E$. There is no need to separately consider $H$, or $M$, once that is done, while $m$ only arises as the magnitude of the 4-vector formed from $E$ and $𝐩$: $\left(mc^2\right)^2 = E^2 - |𝐩c|^2$. So, it doesn't need to be independently considered, either.
Thus, the stress tensor and stress tensor density are $4×4$: in some fashion or another it contains the densities for $\left(E, 𝐩\right)$ and their respective fluxes (each of which is a 3-vector). That adds up to 16 components in all.
Let $(M,L,T)$ denote respectively the dimension of mass, length and time duration. On a coordinate grid patterned after a Cartesian grid, where all spatial coordinates have the dimension $L$, and where the time coordinate has the dimension $T$, the dimensions are the following, where $H = ML^2/T$ the dimension of action:
$$
\left[𝔗^0_0\right] = \frac{H}{L^3T} = \frac{M}{LT^2} = \frac{1}{L^3}\frac{ML^2}{T^2}, \hspace 1em \left[𝔗^i_0\right] = \left[𝔗^0_0\right]\frac{L}{T} \hspace 1em (i = 1, 2, 3), \\
\left[𝔗^0_j\right] = \frac{H}{L^3T}\frac{T}{L} = \frac{M}{L^2T} = \frac{1}{L^3}\frac{ML}{T}, \hspace 1em \left[𝔗^i_j\right] = \left[𝔗^0_j\right]\frac{L}{T} \hspace 1em (i,j = 1, 2, 3).
$$
In order, from top to bottom, left to right, that's: energy density and its flux, for $p_0 = -E$; momentum density and its flux, for $\left(p_1, p_2, p_3\right) = 𝐩$. The dimension for momentum flux is the same as that for pressure and stress, since:
$$\frac{1}{L^3}\frac{ML}{T}\frac{L}{T} = \frac{M}{LT^2} = \frac{1}{L^2}\frac{ML}{T^2}.$$
It's also the same as the dimension for energy density.
Technically, if you're doing this on grids based on cylindrical or spherical coordinates, where some of the spatial coordinates are angular (and dimensionless), then that changes things, since you're actually supposed to include the Jacobian in the density, if you're referencing to angular coordinates. For a cylinder with radius $r$, that's $r$, while for spherical coordinates with co-latitude $θ$ (which is 0 degrees at the north pole, 180 degrees at the south pole), it's $r^2\sin{θ}$.
In Newtonian theory, what was 4 remains as 5: the mass $M = m$, momentum $𝐩$ and kinetic energy $H$. There is no Newtonian analogue of the mass-energy unification of $H$ and $M$ into $E$: no reduction $\left(M, 𝐩, H\right)$ into $\left(E, 𝐩\right)$.
So, what does its stress tensor (density) look like?! And what corresponds to what? Things don't exactly line up and need to be fixed.
Actually, the fix needed here is with Relativity.
What became 4, in Relativity, can still be kept as 5. For a body of mass $m$ and velocity $𝐯$ and speed $v = |𝐯|$, the respective quantities are:
$$𝐩 = M𝐯, \hspace 1em E = Mc^2, \hspace 1em H = Mw.$$
What is $w$? Well, that comes out of the equations expressing the relation between these quantities:
$$M^2 = m^2 + \frac{p^2}{c^2}, \hspace 1em\text{sgn}(M) = \text{sgn}(m), \hspace 1em M = m + \frac{H}{c^2}.$$
What we find is that
$$
M = \frac{m}{\sqrt{1 - (v/c)^2}},\\
H = c^2(M - m) = mc^2\left(\frac{1}{\sqrt{1 - (v/c)^2}} - 1\right) =
\frac{m}{{\sqrt{1 - (v/c)^2}}}\frac{v^2}{1 + \sqrt{1 - (v/c)^2}},
$$
so $w$ is given by:
$$w = \frac{v^2}{1 + \sqrt{1 - (v/c)^2}} = c^2\left(1 - \sqrt{1 - (v/c)^2}\right).$$
In Relativity, that's the Lagrangian per unit mass for a body, after removing the $c^2$ offset.
The equation $E^2 - |𝐩c|^2 = \left(mc^2\right)^2$ now splits into two equations:
$$|𝐩|^2 - 2MH + \frac{1}{c^2}H^2 = 0, \hspace 1em M - \frac{1}{c^2}H = m. \tag{Rel1}$$
Correspondingly, we also have:
$$v^2 - 2w + \frac{1}{c^2}w^2 = 0. \label{Rel2}\tag{Rel2}$$
This carries with it, the strong suggestion to treat $w$, itself, as a velocity component! We'll write it as:
$$w = -\frac{du}{dt}.$$
To what end? This: all the above has a non-relativistic version obtained by simply setting $(1/c)^2 = 0$:
$$|𝐩|^2 - 2MH = 0, \hspace 1em M = m, \tag{New1}$$
and
$$v^2 - 2w = 0. \label{New2}\tag{New2}$$
If we also write the velocity as
$$𝐯 = \frac{d𝐫}{dt}, \hspace 1em 𝐫 = (x, y, z),$$
then from the velocity relations ($\ref{Rel2}$) and ($\ref{New2}$), we can write down the corresponding geometric invariants which also double over as constraints:
$$|d𝐫|^2 + 2 dt du + \frac{1}{c^2}du^2 = 0, \label{Rel3}\tag{Rel3}$$
and
$$|d𝐫|^2 + 2 dt du = 0. \tag{New3}$$
These are the flat, zero-gravity, limits of geometries respectively corresponding to Relativity and Newtonian physics. The latter is called Bargmann geometry. The former is the (unnamed) Relativistic version of it.
What is $u$ and what relation to Minkowski geometry does ($\ref{Rel3}$) have? The answer is contained in the identity:
$$1 - \frac{w}{c^2} = \sqrt{1 - \left(\frac{v}{c}\right)^2},$$
which allows us to an equation for the proper time, which we'll denote $ds$:
$$ds ≡ dt \sqrt{1 - \left(\frac{v}{c}\right)^2} = dt \left(1 - \frac{w}{c^2}\right) = dt + \frac{1}{c^2} du. \tag{Rel4}$$
The Newtonian analogue of this - again obtained by setting $(1/c)^2 = 0$ - merely states that proper time is ordinary time:
$$ds = dt. \tag{New4}$$
In Relativity, $u = c^2(s - t)$ is just the time dilation, itself, ramped up by a factor of $c^2$. It continues to live on in the non-relativistic setting, though its moorings to "time dilation" have been cut loose. Its velocity $du/dt = -w = -v^2/2$ involves the same $w = v^2/2$ that appears in the formula for kinetic energy: $H = Mw$.
A more direct answer can be arrived at by noting that if we do, indeed, treat $H$ as the Hamiltonian for the body, then the corresponding Lagrangian would be expressed as:
$$L = 𝐩·𝐯 - H = M(v^2 - w) = mw.$$
Thus, $w$ is the Lagrangian per unit mass for the body, while $u = -\int L/m dt = -S/m$ is the action per unit mass, apart from the sign.
Substituting $du$ for $ds$ back into ($\ref{Rel3}$), one obtains the following equation for the proper time:
$$ds^2 = dt^2 - \frac{|d𝐫|^2}{c^2}. \tag{Rel5}$$
There is no Newtonian analogue of this - no reduction of the 5D Bargmann geometry to any 4D geometry that has a self-contained, invertible, metric and line element.
So, now ... what happens if you make the substitution for $ds$ by $du$ in the Schwarzschild solution? First, rewrite it as a proper time metric:
$$-c^2ds^2 = \frac{dr^2}{1 - k/r} + r^2\left(dθ^2 + (\sin θ dφ)^2\right) - c^2 dt^2 \left(1 - \frac{k}{r}\right),$$
and bring everything over to one side of the equation
$$c^2ds^2 + dx^2 + dy^2 + dz^2 + \left(\frac{1}{1 - k/r} - 1\right)dr^2 - c^2dt^2 + \frac{c^2k}{r}dt^2 = 0,$$
after rewriting the part involving the angular coordinates back into Cartesian form:
$$r^2\left(dθ^2 + (\sin θ dφ)^2\right) = dx^2 + dy^2 + dz^2 - dr^2.$$
After substituting in $du$ and making further simplifications, we get:
$$dx^2 + dy^2 + dz^2 + 2 dt du + \frac{1}{c^2}du^2 + \frac{c^2k}{r}dt^2 + \left(\frac{k}{r - k}\right)dr^2 = 0.$$
There is no meaningful non-relativistic version of this, nor of ($\ref{Rel0}$), unless we assume that the model parameter $k$ scales as $(1/c)^2$. So, we'll rewrite it as $k = K/c^2$ and rewrite the line element as:
$$dx^2 + dy^2 + dz^2 + 2 dt du + \frac{1}{c^2}du^2 + \frac{K}{r}dt^2 + \frac{1}{c^2}\left(\frac{K}{r - K/c^2}\right)dr^2 = 0, \\ ds = dt + \frac{1}{c^2}du. \tag{Rel0A}$$
The Newtonian version of this is:
$$dx^2 + dy^2 + dz^2 + 2 dt du + \frac{K}{r}dt^2 = 0, \hspace 1em ds = dt. \tag{New0A}$$
The corresponding constraint for velocity is:
$$v^2 - 2w + \frac{K}{r} = 0 \hspace 1em⇒\hspace1em w = \frac{v^2}{2} + \frac{K}{2r}.$$
That's not an approximation or anything. That's for real. The geodesics on this geometry will give the orbital motion for Newton's law of gravity, for the right setting of $K$; plus $u$ will come out to be negative of the total action per unit mass for the orbiting body. This is the line element for Newtonian Gravity - as a warped Bargmann geometry. It is a 4+1 dimensional geometry, and geodesic motions on it are constrained to be on its light cone. Everything Newtonian is "light speed" in 5D Bargmann space.
For a body orbiting a big spherical clump of mass, like a planet, that has mass $μ$, the body's Lagrangian per unit mass will be given by:
$$\frac{L}{m} = \frac{v^2}{2} + \frac{Gμ}{r},$$
which establishes the values of the model parameters:
$$K = 2Gμ, \hspace 1em k = \frac{2Gμ}{c^2}.$$
So, now we get to the matter of the interior solution. Assume the sphere is of radius $R > 0$ and that the matter density $ρ$ is uniform throughout. That means that the total mass will be given by:
$$μ = \frac{4πρR^3}{3}.$$
For the interior solution, we will continue to use ($\ref{Rel0}$) for $0 < r < R$, by treating the mass shell that lies in the range $r ≤ R$ as a "Faraday Cage" for gravity, meaning: its effect cancels out. The corresponding mass will then be $μ_r = μ(r/R)^3$. Likewise, $k$ will now be a function of $r$ as follows:
$$k(r) = \frac{2Gμ_r}{c^2} = \frac{2Gμr^3}{R^3c^3}, \hspace 1em ρ = \frac{3μ}{4πR^3}.$$
The line element in ($\ref{Rel0}$) can then be rewritten as:
$$\frac{dr^2}{A} + r^2\left(dθ^2 + (\sin θ dφ)^2\right) - c^2 A dt^2, \hspace 1em A = 1 - \frac{k(r)}{r}. \label{Rel0r}\tag{Rel0r}$$
So, now it comes down to the matter of substituting ($\ref{Rel0r}$) on the left-hand side of Einstein's equation - but written for the stress tensor density: $\sqrt{|g|}g^{ρμ} G_{μν} = κ 𝔗^ρ_ν$ and then noting that $𝔗^0_0$ is the density for the component $p_0 = -E = -Mc^2$ of the momentum 4-vector, we may write $\sqrt{|g|}g^{0μ} G_{μ0} = κ 𝔗^0_0 = -κρc^2$, if working in coordinates patterned on Cartesian grids, where $ρ$ is the density for the moving mass $M$ (not the rest mass $m$). The distinction, here, between $m$ and $M$ is actually moot, since the gravitating body is stationary.
As noted above, if you're not working in coordinates patterned on a Cartesian grid, but in coordinates patterned on cylindrical or spherical grid, then the integral for mass, in terms of density $ρ$ will actually have the form:
$$μ_r = \int_0^r ρ r^2 \sin{θ} dr dθ dφ.$$
So, technically: the density when referred to a spherical grid would be $𝔗^0_0 = -ρc^2 r^2 \sin{θ}$ - and that is what we will actually be matching against! So, we should get an equation of the form:
$$\sqrt{|g|} g^{0μ} G_{μ0} = κ 𝔗^0_0 = -κρc^2 r^2 \sin{θ}.$$.
The dimensional analysis can be done for any type of coordinate grid, the bottom line will be the same. It is easiest for a Cartesian-based grid, where the dimensional analysis runs as follows:
$$\left[G_{00}\right] = \frac{1}{T^2}, \hspace 1em [ρ] = \frac{M}{L^3}, \hspace 1em [c] = \frac{L}{T}, \hspace 1em \left[1 - \frac{k}{r}\right] = 1.$$
So, we must have:
$$[κ] = \frac{1/T^2}{M/T^3} = \frac{T}{M}.$$
Since Poisson's equation - which is the Newtonian limit of this $∇^2φ = 4πGρ$, (where $φ$ is the gravitational potential per unit mass) - only has $ρ$ in combination with Newton's coefficient $G$ then $ρ$ may only appear in the combination $Gρ$, here, which means $κ = \left(\_\right)·G$ where the other factors are independent of $G$. The dimension of $G$ is $[G] = L^3/(MT^2)$. Therefore, the ratio must have dimension:
$$\left[\frac{κ}{G}\right] = \frac{T/M}{L^3/(MT^2)} = \frac{T^3}{L^3} = \left[\frac{1}{c^3}\right].$$
Thus, $κ = \left(\_\right)·G/c^3$, where the remaining factor is numerical and dimensionless.
So, with our "Faraday Cage" heuristic, after substituting in the left hand side and (not) doing a long calculation (because the cheat sheet over on the Wikipedia spares us the labor), we get the following for the Ricci tensor (after substituting $A = 1 - k(r)/r$ in for $R_{22}$ and $R_{33}$):
$$R_{00} = \frac{c^2}{2}A∇^2A, \hspace 1em R_{11} = -\frac{1}{2}\frac{∇^2A}{A}, \hspace 2em R_{22} = 1 - \frac{d}{dr}(rA) = \frac{dk}{dr}, \hspace 1em R_{33} = R_{22} \sin^2{θ},$$
all other components being zero. After using the metric to raise the indices, we get a diagonal matrix with the following components:
$$R^0_0 = -\frac{∇^2A}{2} = R^1_1, \hspace 1em R^2_2 = \frac{dk/dr}{r^2}$$
After subtracting out half the trace
$$R = R^ρ_ρ = -∇^2A + 2\frac{dk/dr}{r^2},$$
we get the Einstein tensor $G^ρ_ν = R^ρ_ν - ½ δ^ρ_ν R$:
$$G^0_0 = -\frac{dk/dr}{r^2} = G^1_1, \hspace 1em G^2_2 = \frac{∇^2A}{2} = G^3_3.$$
Finally - and most importantly - is the distinction between tensors and tensor densities! We're actually dealing with densities. So, after densitizing, i.e. multiplying by $\sqrt{|g|} = c r^2 \sin{θ}$, and writing $𝔊^ρ_ν = \sqrt{|g|} G^ρ_ν$, we get:
$$𝔊^0_0 = -c \frac{dk}{dr} \sin{θ} = 𝔊^1_1, \hspace 1em 𝔊^2_2 = \frac{c}{2}∇^2Ar^2\sin{θ} = 𝔊^3_3.$$
It is the component $𝔊^0_0$ we are interested in. Substitute in for $k(r)$ to get:
$$\frac{dk}{dr} = \frac{8πGρr^2}{c^2}.$$
Thus,
$$𝔊^0_0 = -\frac{8πGρ}{c} r^2\sin{θ}.$$
Equate this to $κ𝔗^0_0 = -κρc^2r^2\sin{θ}$, we get the following match:
$$𝔊^0_0 = κ𝔗^0_0 \hspace 1em⇒\hspace 1em -\frac{8πGρ}{c}r^2\sin{θ} = -κρc^2r^2\sin{θ} \hspace 1em⇒\hspace 1em κ = \frac{8πG}{c^3}.$$
Thus, the coefficient is fixed.
I left a loose end dangling: is there a Newtonian analogue to Einstein's equation? So, then is the stress tensor density $5×5$? Is it a symmetric tensor, when it's de-densitized and its upper index is lowered, using the (Bargmann) metric? The answer is yes to all of the above. It is not what Cartan wrote in the 1920's, but much cleaner. But ... it doesn't fully integrate fluid dynamics into the law of gravity (neither did Cartan's version). That gets lost in the non-relativistic limit $(1/c)^2 → 0$, because there are residual terms of the order $(1/c)^4$ whose survival into the non-relativistic limit have to be ensured, in order to get the fluid-dynamics / gravity unification that Einstein accomplishes with his equation.
Appendix: Some Details On The Calculation
For the metric:
$$
g_{00} = -c^2A, \hspace 1em g_{11} = \frac{1}{A}, \hspace 1em g_{22} = r^2, \hspace 1em g_{33} = (r \sin{θ})^2, \hspace 1em g_{μν} = 0 \hspace 1em (μ ≠ ν), \\
A = 1 - \frac{k(r)}{r}, \hspace 1em g = -\left(c r^2 \sin{θ}\right)^2, \hspace 1em \sqrt{|g|} = c r^2 \sin{θ}
$$
For the connection:
$$
Γ^1_{00} = +½c^2A'A, \hspace 1em Γ^0_{01} = Γ^0_{10} = +½\frac{A'}{A} \\
Γ^1_{11} = -½\frac{A'}{A}, \\
Γ^2_{12} = Γ^2_{21} = \frac{1}{r}, \hspace 1em Γ^1_{22} = -rA, \\
Γ^3_{13} = Γ^3_{31} = \frac{1}{r}, \hspace 1em Γ^1_{33} = -rA \sin^2{θ}, \\
Γ^3_{23} = Γ^3_{32} = \cot{θ}, \hspace 1em Γ^2_{33} = -\sin{θ}\cos{θ}
$$
For the Ricci coefficients, scalar coefficient and Einstein tensor:
$$
R_{00} = ½ c^2A ∇^2A, \hspace 1em R_{11} = -½\frac{∇²A}{A}, \hspace 1em R_{22} = k', \hspace 1em R_{33} = k' \sin^2{θ}, \\
R^0_0 = -½ ∇^2A = R^1_1, \hspace 1em R^2_2 = \frac{k'}{r^2} = R^3_3, \\
R = -∇^2A + \frac{2k'}{r^2}, \\
G^0_0 = -\frac{k'}{r^2} = G^1_1, \hspace 1em G^2_2 = ½ ∇^2A = G^3_3.
$$
For $k$ and $k'$:
$$k = \frac{8πG}{3} \frac{ρr^3}{c^2} = \frac{8πGρ}{3c^2} r^3 \hspace 1em⇒\hspace 1em k' = \frac{8πG}{c^2} r^2.$$
The densitized Einstein tensor $𝔊^ρ_ν = \sqrt{|g|} G^ρ_ν$:
$$𝔊^0_0 = -\frac{8πGρ}{c} \left(r^2 \sin{θ}\right) = 𝔊^1_1, \hspace 1em 𝔊^2_2 = ½ c ∇^2A \left(r^2 \sin{θ}\right) = 𝔊^3_3.$$