Calculating the $\mu\rightarrow e\gamma$ decay and cancellation between diagrams

(i) It is not advisable to do loop calculations in unitary gauge due to blowing up of the propagator for W boson in ultraviolet limit as the loop involves integrating over momenta. Renormalizability(perturbative) is also obscured by using unitary gauge. So I doubt you will find any good reference on doing loop calculation using unitary gauge in this case.

(ii) A term like $\bar u\gamma_\mu u$ can be extracted from the loop diagrams (a)-(d) by noting the vertex where fermion line is present and fermion propagator contribute one gamma factor. There are three gammas coming with some of them multiplying with $\gamma_5$. some duality relations like $\epsilon_{abcd}\gamma^{bcd}=i\gamma_5\gamma_a$ and simple gamma technology will produce $\bar u\gamma_\mu u$ and $\bar u\gamma_{\mu\nu} u$ type term. [$\gamma_{\mu\nu}$ is same as your $\sigma_{\mu\nu}$(apart from some numerical factor probably)]

(iii) The term $\gamma\centerdot\epsilon$ does not cancel from the diagrams but is rejected on the reality condition of emitted photon. Simply saying the term violates gauge invariance. To see it explicitly, $\epsilon^\nu \bar u_e(1+\gamma_5)\gamma_\nu u_\mu$ with $\epsilon^\nu$ removed, operate with $q^\nu$ on $\bar u_e (1+\gamma_5)\gamma_\nu u_\mu$ and write q as q-p+p and then use Dirac equation in momentum space for muon and electron and you will see it does not vanish.

For the edit part- In the massless limit of electron, you can pick a certain helicity and the coupling of W boson is to left handed electron. With the bar form of electron spinor, it basically requires to multiply it with $1+\gamma_5$. This amounts to $(1+\gamma_5)(A+B\gamma_5)$ or $A(1+\gamma_5)+B(1+\gamma_5)$ as $\gamma_5$ commutes with $\sigma_{\mu\nu}$. So both term involving A and B corresponds to same matrix element in massless electron limit, So they can be taken to be equal .