[proofplan]
We prove the multiplier identity by computing the distributional Fourier transform of the principal-value kernel $K_j(x) = c_n x_j |x|^{-(n+1)}$. The kernel $K_j$ is odd and homogeneous of degree $-n$, so by general principles its Fourier transform is even and homogeneous of degree $0$ — an even bounded function on the unit sphere. We identify the explicit form by exploiting the relation between $K_j$ and the partial derivative of the Riesz potential $|x|^{1-n}$: writing $K_j = -c_n'\, \partial_{x_j}|x|^{1-n}$ for an explicit constant, and using the known Fourier transform of $|x|^{1-n}$, we obtain $\hat{K}_j(\xi) = -i\xi_j/|\xi|$. The constant $c_n = \Gamma((n+1)/2)/\pi^{(n+1)/2}$ is normalised exactly so that the prefactor cleans up to $-i\xi_j/|\xi|$.
[/proofplan]
[step:Express $K_j$ as a distributional derivative of the Riesz potential]
Define the Riesz potential function
\begin{align*}
I: \mathbb{R}^n \setminus \{0\} &\to \mathbb{R} \\
x &\mapsto |x|^{1-n}.
\end{align*}
Since $I$ is locally integrable on $\mathbb{R}^n$ (the singularity at the origin has integrable size: $\int_{B(0,1)} |x|^{1-n}\, d\mathcal{L}^n(x) = n\alpha_n \int_0^1 r^{1-n}\cdot r^{n-1}\, d\mathcal{L}^1(r) = n\alpha_n < \infty$, where $\alpha_n$ is the volume of the unit ball), $I \in L^1_{\mathrm{loc}}(\mathbb{R}^n)$, and $I$ defines a tempered distribution via
\begin{align*}
T_I(\varphi) := \int_{\mathbb{R}^n} I(x)\, \varphi(x)\, d\mathcal{L}^n(x), \qquad \varphi \in \mathcal{S}(\mathbb{R}^n).
\end{align*}
For $x \neq 0$, the classical partial derivative is
\begin{align*}
\partial_{x_j} |x|^{1-n} = (1-n)|x|^{-n-1}x_j = -(n-1)\frac{x_j}{|x|^{n+1}}.
\end{align*}
We claim that the distributional partial derivative of $T_I$ equals the principal-value distribution
\begin{align*}
P_j(\varphi) := -(n-1)\,\mathrm{p.v.}\!\int_{\mathbb{R}^n} \frac{x_j}{|x|^{n+1}}\varphi(x)\, d\mathcal{L}^n(x).
\end{align*}
For $\varepsilon > 0$, set $I_\varepsilon(x) := |x|^{1-n}\mathbb{1}_{\{|x|>\varepsilon\}}$. Then $I_\varepsilon \in L^1_{\mathrm{loc}}(\mathbb{R}^n)$ has a continuous classical partial derivative on $\{|x| > \varepsilon\}$, and integrating by parts over $\{|x| > \varepsilon\}$ against $\varphi \in \mathcal{S}(\mathbb{R}^n)$,
\begin{align*}
-\int_{|x|>\varepsilon} I(x)\, \partial_{x_j}\varphi(x)\, d\mathcal{L}^n(x) = -(n-1)\int_{|x|>\varepsilon} \frac{x_j}{|x|^{n+1}}\varphi(x)\, d\mathcal{L}^n(x) - \int_{|x|=\varepsilon} I(x)\,\varphi(x)\,\nu_j(x)\, d\mathcal{H}^{n-1}(x),
\end{align*}
where $\nu(x) = -x/|x|$ is the outward normal to the region $\{|x|>\varepsilon\}$ on the sphere $\{|x|=\varepsilon\}$, so $\nu_j(x) = -x_j/|x| = -x_j/\varepsilon$. The boundary term equals
\begin{align*}
-\int_{|x|=\varepsilon} \varepsilon^{1-n}\,\varphi(x)\cdot\left(-\frac{x_j}{\varepsilon}\right)\, d\mathcal{H}^{n-1}(x) = \varepsilon^{-n}\int_{|x|=\varepsilon} x_j\,\varphi(x)\, d\mathcal{H}^{n-1}(x).
\end{align*}
As $\varepsilon \to 0^+$, $\int_{|x|=\varepsilon} x_j\, \varphi(x)\, d\mathcal{H}^{n-1}(x) = \varepsilon^n \int_{S^{n-1}} \omega_j\, \varphi(\varepsilon\omega)\, d\sigma(\omega)$, where $\sigma$ is the surface measure on the unit sphere $S^{n-1}$. The integrand $\omega_j \varphi(\varepsilon\omega) \to \omega_j\varphi(0)$ pointwise, and $\int_{S^{n-1}}\omega_j\, d\sigma(\omega) = 0$ by symmetry. So the boundary term vanishes as $\varepsilon \to 0^+$:
\begin{align*}
\varepsilon^{-n}\int_{|x|=\varepsilon} x_j \varphi(x)\, d\mathcal{H}^{n-1}(x) \to \int_{S^{n-1}}\omega_j\, d\sigma(\omega) \cdot \varphi(0) = 0.
\end{align*}
On the other hand, the LHS converges to $-T_I(\partial_{x_j}\varphi) = (\partial_{x_j} T_I)(\varphi)$ by dominated convergence (with dominator $|I|\cdot |\partial_{x_j}\varphi|$, which is in $L^1$ since $\partial_{x_j}\varphi \in \mathcal{S}$). The RHS converges to $P_j(\varphi)$ by the definition of the principal value. Hence
\begin{align*}
\partial_{x_j} T_I = P_j \quad \text{in } \mathcal{S}'(\mathbb{R}^n).
\end{align*}
In terms of the kernel $K_j(x) := c_n\, x_j/|x|^{n+1}$, we have $P_j = -\frac{n-1}{c_n}\cdot c_n\, \mathrm{p.v.}\!\int K_j \varphi$, so the principal-value tempered distribution associated with $K_j$ — call it $T_{K_j}$ — satisfies
\begin{align*}
T_{K_j} = -\frac{1}{n-1}\cdot \partial_{x_j}T_{c_n I}, \qquad T_{c_n I}(\varphi) = c_n \int_{\mathbb{R}^n} |x|^{1-n}\varphi(x)\, d\mathcal{L}^n(x).
\end{align*}
Equivalently, $T_{K_j} = -\frac{c_n}{n-1}\, \partial_{x_j}T_I$.
[/step]
[step:Recall the Fourier transform of the Riesz potential]
We invoke the following known identity for the Fourier transform of the Riesz potential. Let $0 < \alpha < n$. Then $|x|^{-\alpha} \in L^1_{\mathrm{loc}}(\mathbb{R}^n)$ defines a tempered distribution, and
\begin{align*}
\mathcal{F}(|x|^{-\alpha})(\xi) = \gamma_{n,\alpha}\, |\xi|^{\alpha - n}, \qquad \gamma_{n,\alpha} = \frac{2^{n-\alpha}\pi^{n/2}\Gamma((n-\alpha)/2)}{\Gamma(\alpha/2)},
\end{align*}
where the equality is in $\mathcal{S}'(\mathbb{R}^n)$. See [Fourier Transform of Riesz Potentials](/theorems/???) for the derivation; the proof goes via analytic continuation or via the Gaussian heat kernel representation $|x|^{-\alpha} = \frac{1}{\Gamma(\alpha/2)}\int_0^\infty t^{\alpha/2-1} e^{-t|x|^2}\, d\mathcal{L}^1(t)$ (after rescaling).
Applying this with $\alpha = n - 1$, so that $|x|^{-\alpha} = |x|^{1-n}$ and $|\xi|^{\alpha - n} = |\xi|^{-1}$, gives
\begin{align*}
\mathcal{F}(|x|^{1-n})(\xi) = \gamma_{n,n-1}\, |\xi|^{-1}, \qquad \gamma_{n, n-1} = \frac{2 \pi^{n/2}\,\Gamma(1/2)}{\Gamma((n-1)/2)} = \frac{2\pi^{(n+1)/2}}{\Gamma((n-1)/2)},
\end{align*}
using $\Gamma(1/2) = \sqrt{\pi}$.
[/step]
[step:Differentiate on the Fourier side and combine constants]
Distributional differentiation on the spatial side becomes multiplication by $i\xi_j$ on the Fourier side:
\begin{align*}
\mathcal{F}(\partial_{x_j} T)(\xi) = i\xi_j\, \mathcal{F}(T)(\xi), \qquad T \in \mathcal{S}'(\mathbb{R}^n).
\end{align*}
Applying this to $T = T_{|x|^{1-n}}$,
\begin{align*}
\mathcal{F}(\partial_{x_j}T_{|x|^{1-n}})(\xi) = i\xi_j \cdot \frac{2\pi^{(n+1)/2}}{\Gamma((n-1)/2)|\xi|} = \frac{2\pi^{(n+1)/2}}{\Gamma((n-1)/2)}\cdot \frac{i\xi_j}{|\xi|}.
\end{align*}
By Step 1, $T_{K_j} = -\frac{c_n}{n-1}\partial_{x_j}T_{|x|^{1-n}}$. Therefore
\begin{align*}
\hat{K}_j(\xi) = \mathcal{F}(T_{K_j})(\xi) = -\frac{c_n}{n-1}\cdot \frac{2\pi^{(n+1)/2}}{\Gamma((n-1)/2)}\cdot\frac{i\xi_j}{|\xi|}.
\end{align*}
We now compute the constant. Recall $c_n = \Gamma((n+1)/2)/\pi^{(n+1)/2}$. The reflection identity for the Gamma function $\Gamma((n+1)/2) = \frac{n-1}{2}\Gamma((n-1)/2)$ gives
\begin{align*}
\frac{c_n}{n-1}\cdot \frac{2\pi^{(n+1)/2}}{\Gamma((n-1)/2)} = \frac{\Gamma((n+1)/2)}{\pi^{(n+1)/2}(n-1)}\cdot \frac{2\pi^{(n+1)/2}}{\Gamma((n-1)/2)} = \frac{2\Gamma((n+1)/2)}{(n-1)\Gamma((n-1)/2)} = \frac{2 \cdot \frac{n-1}{2}\Gamma((n-1)/2)}{(n-1)\Gamma((n-1)/2)} = 1.
\end{align*}
Hence
\begin{align*}
\hat{K}_j(\xi) = -\frac{i\xi_j}{|\xi|}.
\end{align*}
[/step]
[step:Conclude the multiplier identity for the Riesz transform]
For $f \in \mathcal{S}(\mathbb{R}^n)$, the principal-value definition of the Riesz transform reads $R_j f = K_j * f$, the convolution of the tempered distribution $T_{K_j}$ with the Schwartz function $f$. By the convolution theorem for the Fourier transform of $\mathcal{S}' \times \mathcal{S}$,
\begin{align*}
\widehat{R_j f}(\xi) = \widehat{K_j * f}(\xi) = \hat{K}_j(\xi) \cdot \hat{f}(\xi) = -i\frac{\xi_j}{|\xi|}\hat{f}(\xi),
\end{align*}
which is the claimed identity. The Fourier transform of $K_j * f$ is well-defined because the multiplier $-i\xi_j/|\xi|$ is bounded by $1$ on $\mathbb{R}^n \setminus \{0\}$, so multiplication by it preserves $L^2$ and the resulting tempered distribution coincides with the function $R_j f \in L^2(\mathbb{R}^n)$.
[/step]