Proof of Cauchy Riemann Equations in Polar Coordinates
I happen to have some notes on this question. What follows here is the usual approach, it's just multivariate calculus paired with the Cauchy Riemann equations. I have an idea for an easier way, I'll post it as a second answer in a bit if it works.
If we use polar coordinates to rewrite $f$ as follows: $$ f(x(r,\theta),y(r,\theta)) = u(x(r,\theta),y(r,\theta))+iv(x(r,\theta),y(r,\theta)) $$ we use shorthands $F(r,\theta)=f(x(r,\theta),y(r,\theta))$ and $U(r,\theta )=u(x(r,\theta),y(r,\theta))$ and $V(r,\theta )=v(x(r,\theta),y(r,\theta))$. We derive the CR-equations in polar coordinates via the chain rule from multivariate calculus, $$ U_r = x_ru_x + y_ru_y = \cos(\theta)u_x + \sin(\theta)u_y \ \ \text{and} \ \ U_{\theta} = x_{\theta}u_x + y_{\theta}u_y = -r\sin(\theta)u_x + r\cos(\theta)u_y $$ Likewise, $$ V_r = x_rv_x + y_rv_y = \cos(\theta)v_x + \sin(\theta)v_y \ \ \text{and} \ \ V_{\theta} = x_{\theta}v_x + y_{\theta}v_y = -r\sin(\theta)v_x + r\cos(\theta)v_y $$ We can write these in matrix notation as follows: $$ \left[ \begin{array}{l} U_r \\ U_{\theta} \end{array} \right] = \left[ \begin{array}{ll} \cos(\theta) & \sin(\theta) \\ -r\sin(\theta) & r\cos(\theta) \end{array} \right]\left[ \begin{array}{l} u_x \\ u_y \end{array} \right] \ \ \text{and} \ \ \left[ \begin{array}{l} V_r \\ V_{\theta} \end{array} \right] = \left[ \begin{array}{ll} \cos(\theta) & \sin(\theta) \\ -r\sin(\theta) & r\cos(\theta) \end{array} \right]\left[ \begin{array}{l} v_x \\ v_y \end{array} \right] $$ Multiply these by the inverse matrix: $\left[ \begin{array}{ll} \cos(\theta) & \sin(\theta) \\ -r\sin(\theta) & r\cos(\theta) \end{array} \right]^{-1} = \frac{1}{r}\left[ \begin{array}{ll} r\cos(\theta) & -\sin(\theta) \\ r\sin(\theta) & \cos(\theta) \end{array} \right]$ to find $$ \left[ \begin{array}{l} u_x \\ u_y \end{array} \right] = \frac{1}{r}\left[ \begin{array}{ll} r\cos(\theta) & -\sin(\theta) \\ r\sin(\theta) & \cos(\theta) \end{array} \right]\left[ \begin{array}{l} U_r \\ U_{\theta} \end{array} \right] = \left[ \begin{array}{l} \cos(\theta)U_r - \tfrac{1}{r}\sin(\theta)U_{\theta} \\ \sin(\theta)U_r + \tfrac{1}{r}\cos(\theta)U_{\theta} \end{array} \right] $$ A similar calculation holds for $V$. To summarize: $$ u_x = \cos(\theta)U_r - \tfrac{1}{r}\sin(\theta)U_{\theta} \ \ \ \ v_x = \cos(\theta)V_r - \tfrac{1}{r}\sin(\theta)V_{\theta} $$ $$ u_y =\sin(\theta)U_r + \tfrac{1}{r}\cos(\theta)U_{\theta} \ \ \ \ v_y =\sin(\theta)V_r + \tfrac{1}{r}\cos(\theta)V_{\theta} $$ The CR-equation $u_x=v_y$ yields: $$ (A.) \ \ \cos(\theta)U_r - \tfrac{1}{r}\sin(\theta)U_{\theta} = \sin(\theta)V_r + \tfrac{1}{r}\cos(\theta)V_{\theta} $$ Likewise the CR-equation $u_y=-v_x$ yields: $$ (B.) \ \ \sin(\theta)U_r + \tfrac{1}{r}\cos(\theta)U_{\theta} = -\cos(\theta)V_r + \tfrac{1}{r}\sin(\theta)V_{\theta}$$ Multiply (A.) by $r\sin(\theta)$ and $(B.)$ by $r\cos(\theta)$ and subtract (A.) from (B.): $$ \boxed{U_{\theta} = -rV_r} $$ Likewise multiply (A.) by $r\cos(\theta)$ and $(B.)$ by $r\sin(\theta)$ and add (A.) and (B.): $$ \boxed{rU_r = V_{\theta}} $$ Finally, recall that $z = re^{i\theta}=r(\cos(\theta)+i\sin(\theta))$ hence \begin{align} \notag f'(z) &= u_x+iv_x \\ \notag &= (\cos(\theta)U_r - \tfrac{1}{r}\sin(\theta)U_{\theta})+i(\cos(\theta)V_r - \tfrac{1}{r}\sin(\theta)V_{\theta}) \\ \notag &= (\cos(\theta)U_r + \sin(\theta)V_{r})+i(\cos(\theta)V_r - \sin(\theta)U_{r}) \\ \notag &= (\cos(\theta)- i\sin(\theta))U_r + i(\cos(\theta)-i\sin(\theta))V_r \\ \notag &= e^{-i\theta}( U_r+iV_r) \notag \end{align}
We can derive using purely polar coordinates. Start with \begin{align} z(r,\theta) &=r\,\mathrm{e}^{\mathrm{i}\theta} \\ f(z) &= u(r,\theta) + \mathrm{i} v(r,\theta) \end{align} We define $f'(z)$ using the limit $$ f'(z) = \lim_{z\to 0} \frac{\Delta f}{\Delta z} $$ where \begin{align} \Delta f &= \Delta u + \mathrm{i} \Delta v \\ \Delta z &= z(r+\Delta r,\theta + \Delta\theta) - z(r,\theta) \\ &= (r+\Delta r)\,\mathrm{e}^{\mathrm{i}(\theta+\Delta\theta)} - r \,\mathrm{e}^{\mathrm{i}\theta} \end{align} Next, we try to first approach from $\Delta\theta\to 0$, $$ \Delta z = (r+\Delta r)\,\mathrm{e}^{\mathrm{i}\theta} - r \,\mathrm{e}^{\mathrm{i}\theta} = \Delta r\,\mathrm{e}^{\mathrm{i}\theta} $$ Therefore when we take the limit, $$ f'(z) = \lim_{\Delta r\to 0} \frac{\Delta u + \mathrm{i} \Delta v}{\Delta r \,\mathrm{e}^{i\theta}} = \frac{1}{\mathrm{e}^{\mathrm{i}\theta}}(u_r + \mathrm{i} v_r) \label{a}\tag{1}$$ On the other hand, approaching from $\Delta r\to 0$ first yields $$\Delta z = r\,\mathrm{e}^{\mathrm{i}(\theta+\Delta\theta)} - r \,\mathrm{e}^{\mathrm{i}\theta}$$ Using the derivative $$\mathrm{e}^{\mathrm{i}(\theta+\Delta\theta)}-\mathrm{e}^{\mathrm{i}\theta} = \frac{\mathrm{d}\mathrm{e}^{\mathrm{i}\theta}}{\mathrm{d}\theta}\Delta\theta = \mathrm{i}\mathrm{e}^{\mathrm{i}\theta}\Delta\theta$$ we get $$\Delta z = \mathrm{i}r\,\mathrm{e}^{\mathrm{i}\theta}\Delta\theta$$ Therefore $$f'(z) = \lim_{\Delta\theta\to 0} \frac{\Delta u + \mathrm{i} \Delta v}{\mathrm{i}r\,\mathrm{e}^{\mathrm{i}\theta}\Delta\theta} = \frac{1}{\mathrm{i}r\,\mathrm{e}^{\mathrm{i}\theta}}(u_\theta + \mathrm{i}v_\theta) \label{b}\tag{2}$$ Finally, comparing the real and imaginary parts of ($\ref{a}$) and ($\ref{b}$) gives us what we want: $$u_r=\frac{1}{r}v_\theta \quad,\quad v_r = -\frac{1}{r}u_\theta$$
Reference: Kwong-tin Tang, Mathematical Methods for Engineers and Scientists 1 (Springer, 2007)
A less standard approach: Take $\{e_{r},e_{\theta}\}$ as a basis in polar coordinates then the Jacobian matrix for any function on $\mathbb{R}^2$ has the form: $$ [df] = \left[ \begin{array}{cc} U_r & U_{\theta} \\ V_r & V_{\theta} \end{array} \right] $$ Define complex-differentiability via complex linearity of the differential; $df(vw)=df(v)w$ at the point in question for all $v,w \in \mathbb{C}$. In particular this gives the beautiful formula: $df(v)=df(1)v$; this means the first column of the Jacobian fixes the second by complex multiplication. It is geometrically clear that $\frac{1}{r}e_{\theta} = ie_r$. Observe that $e_{\theta} = ire_r$ thus $df(e_{\theta}) =df(ire_{r})=irdf(e_r)$. On the other hand, $f=U+iV$ and $$ df(e_r) = U_r+iV_r \qquad \& \qquad df(e_{\theta})=f_{\theta}=U_{\theta}+iV_{\theta} $$ Thus, $U_{\theta}+iV_{\theta}=ir(U_r+iV_r)$ and we derive $$ \boxed{ U_{\theta} = -rV_r \qquad \& \qquad V_{\theta}=rU_r. } $$