subtracting the second one from the first one, we get: $$\sum_i \frac{\rm d}{{\rm d}r}(\psi_{i\lambda}' 
\bar{\psi}_{i\mu} -\bar{\psi}_{i\lambda}'\psi_{i\mu}) = 
-\sum_{ij} \psi_{i\mu}\Delta V_{ij} \bar{\psi}_{j\lambda},
\label{eq:three}$$ as all terms including $V_{ij}$ in the right hand side vanish due to symmetry of $V_{ij}$.

Integrating Eq.(\[eq:three\]) from $r=0$ to the asymptotic region and using the usual Wronskian relationship $W[I_{\ell},O_{\ell}] = -2 {\rm i} k$, we get: $$2 {\rm i} k (\bar{S}_{\mu\lambda} - S_{\lambda\mu})=
 \sum_{ij}\int_0^{\infty}\psi_{i \mu} \Delta V_{ij} \bar{\psi}_{j \lambda} 
 {\rm d} r. 
\label{eq:four}$$ Since the coupled channel equations with a symmetrical potential matrix give a symmetrical $S$-matrix $S_{\lambda\mu} = S_{\mu\lambda} $, the left hand side of Eq.(\[eq:four\]) is equal $$2 {\rm i} k (\bar{S}_{\mu\lambda} - S_{\mu\lambda})= 
  2 {\rm i} k \Delta S_{\mu\lambda}$$ Hence we find, reinstating the $2m/\hbar^2$ factor on the right hand side: $$\Delta S_{\lambda \mu} = \frac{{\rm i} m}{\hbar^2 k} \sum_{ij}\int_0^{\infty}
\psi_{i \lambda} \Delta V_{ij} \bar{\psi}_{j \mu} {\rm d} r. 
\label{eq:delS}$$ This expression is valid for any symmetrical finite perturbations $\Delta V_{ij}$ decreasing at $r\to\infty$ sufficiently rapidly. It is easily extended to any system of coupled-channel equations. For small $\Delta V_{ij}$ we make the Born approximation assumption that $\bar{\psi} \sim \psi$ and get: $$\Delta S_{\lambda \mu} = \frac{{\rm i} m}{\hbar^2 k} \sum_{ij}\int_0^{\infty}
\psi_{i \lambda} 