The Kraus representation of a completely depolarising channel I have an answer, but I want to know how to get from the left hand size to the right hand side:
$$
    \frac{1}{d^2} \sum_{i=1}^{d^2} U_i \rho U_i^{\dagger}
  = \operatorname{Tr}[\rho] \frac{I}{d}.
$$
Here the $U_i$ are $d*d$ orthogonal unitary operators.
(Original image here.)
 A: I'll use here the notation from John Watrous' book. Denote the standard orthonormal basis of the $d$-dimensional space as $e_1,\dots,e_d\in\mathbb{C}^d$, let $\operatorname{L}(\mathbb{C}^d)$ denote the set of $d\times d$ matrices (i.e. linear operators), and let $\operatorname{U}(\mathbb{C}^d)$ denote the set of $d\times d$ unitary matrices. The space $\operatorname{L}(\mathbb{C}^d)$ is itself a $d^2$-dimensional complex vector space with inner product given by 
$$
\langle X,Y\rangle = \operatorname{Tr}(X^*Y)
$$
for all $d\times d$ matrices $X$ and $Y$, where $^*$ denotes the conjugate transpose (i.e., $^\dagger$ in your notation). There is an isomorphism between the spaces of $d\times d$ matrices and vectors in $\mathbb{C}^d\otimes\mathbb{C}^d$ that is defined as 
$$
\operatorname{vec}(uv^*) = u\otimes v
$$
for all $u,v\in\mathbb{C}^d$, and extended linearly to all $d\times d$ matrices. It holds that 
$$
\langle \operatorname{vec}(X),\operatorname{vec}(Y)\rangle = \langle X,Y\rangle = \operatorname{Tr}(X^*Y)
$$
for all $d\times d$ matrices $X$ and $Y$. Finally, the Choi representation of a linear map $\Phi:\operatorname{L}(\mathbb{C}^d)\rightarrow\operatorname{L}(\mathbb{C}^d)$ is the operator $J(\Phi)\in\operatorname{L}(\mathbb{C}^d\otimes\mathbb{C}^d)$ defined as 
$$
J(\Phi)= \sum_{a,b=1}^d \Phi(e_ae_b^*)\otimes e_ae_b^*.
$$
Two such linear mappings are the same if and only if they have the same Choi representation. 
Let $\Phi:\operatorname{L}(\mathbb{C}^d)\rightarrow\operatorname{L}(\mathbb{C}^d)$ be the linear map defined as
$$
\Phi(X)=\frac{\operatorname{Tr}(X)}{d} I.
$$
It is clear that $J(\Phi)=\frac{1}{d}I\otimes I$. This is the map on the right-hand side of the equation in your question.
Now let $\bigl\{U_{a,b}\,:\, a,b\in\{1,\dots,d\}\bigr\}\subset\mathrm{U}(\mathbb{C})$ be an arbitrary collection of orthogonal unitary matrices such that 
$$
\operatorname{Tr}(U_{a,b}^*U_{a',b'})=d\, \delta_{a,a'}\delta_{b,b'}.
$$
Define the operator $V\in\operatorname{L}(\mathbb{C}^d\otimes\mathbb{C}^d)$ as 
$$
V = \frac{1}{d}\sum_{a,b=1}^d e_a\otimes e_b \operatorname{vec}(U_{a,b})^*.
$$
It is evident that $V$ is unitary, since it may be checked that 
$$
VV^* = \sum_{a,b=1}^d e_ae_a^*\otimes e_be_b^* = I\otimes I.
$$
It therefore holds that $I\otimes I = V^*V = \frac{1}{d}\sum_{a,b=1}^d\operatorname{vec}(U_{a,b})\operatorname{vec}(U_{a,b})^*$
Define the linear mapping $\Psi:\operatorname{L}(\mathbb{C}^d)\rightarrow \operatorname{L}(\mathbb{C}^d)$ as $$\Psi(X)=\frac{1}{d^2}\sum_{i=1}^{d^2}U_j XU_j^*$$ for all matrices $X\in\operatorname{L}(\mathbb{C}^d)$.
The Choi representation of $\Psi$ is
\begin{align*}
J(\Psi)
&=\frac{1}{d^2}\sum_{a,b=1}^d \operatorname{vec}(U_{a,b})\operatorname{vec}(U_{a,b})^*\\
& = \frac{1}{d} V^*V\\
& = \frac{1}{d} I \otimes I \\
& = J(\Phi),
\end{align*}
and thus $\Phi=\Psi$.
A: The other answer already provides a full proof passing through Choi's representation.
Here I'll give an alternative, direct proof.
The set of matrices $\{U_i\}_{i=1}^{d^2}$ is, by definition, complete.
This means that any matrix $\rho$ can be decomposed with it as (the $d$ normalization factor here follows from the normalization choice $\newcommand{\Tr}{\operatorname{Tr}}\Tr(U_i^\dagger U_j)=d \delta_{ij}$):
$$\rho = \frac{1}{d}\sum_{i=1}^{d^2} \Tr(U_i^\dagger \rho )U_i.\tag A$$
In particular, such decomposition applied to the "component matrices" (1) $|j\rangle\langle k|$ gives:
$$
|j\rangle\langle k|=\frac{1}{d} \sum_{i=1}^{d^2}(u^*_i)_{jk} U_i.
$$
But of course, because $\langle m|j\rangle\langle k| n\rangle=\delta_{mj}\delta_{kn}$, the above immediately gives
$$
\delta_{mj}\delta_{kn} =
\frac{1}{d}\sum_{i=1}^{d^2} (u_i^*)_{jk} (u_i)_{mn}. \tag B
$$
Strong of (B), we may now see that
$$
\frac{1}{d^2}\sum_{i=1}^{d^2} U_i\rho U_i^\dagger =
\frac{1}{d^2}\sum_{i=1}^{d^2} \sum_{jklm} (u_i)_{jk} \rho_{kl} (u^*_i)_{ml} |j\rangle\!\langle m| =
\frac{1}{d} \Tr(\rho) I.
$$

Appendix 1
Interestingly, (B) can also be seen as arising directly from the completeness of the operators $\newcommand{\bU}{\boldsymbol{U}}U_i$.
To see this, let us denote with $\tilde{\bU}_i$ the vectorization of $U_i$, so that if $U_i=\sum_{jk}(u_i)_{jk} |j\rangle\langle k|$, then
$\tilde{\bU}_i\equiv\sum_{jk}(u_i)_{jk}|j\rangle|k\rangle$.
Is is then easy to see that the orthogonality/orthonormality of the operators $U_i$ in the Hilbert-Schmidt inner product is equivalent to the orthogonality/orthonormality of the vectors $\tilde{\bU}_i$ in the regular Euclidean inner product:
$$\Tr(U_i^\dagger U_j)=\langle \tilde{\bU}_i, \tilde{\bU}_j\rangle.$$
The completeness of the $U_i$ is then stated as the completeness relation for the $\tilde{\bU}_i$:
$$\sum_i \tilde{\bU}_i\tilde{\bU}_i^* = I,\tag C$$
where $\tilde{\bU}_i^*$ is the dual of $\tilde{\bU}_i$.
As can be readily checked, (B) is nothing but (C) after expansion of the indices.
Appendix 2: diagrammatic notation
This identity can also be directly shown via diagrammatic notation, if one is familiar with it.
Equation (C) above can be written as

which is equivalent to

Applying the above to a matrix $\rho$, it follows that

which is the result.

(1) I'm here using the convention, common in quantum mechanics, of writing basis vectors as $|j\rangle$ instead of $\boldsymbol{e_j}$, and their duals as $\langle j|$ instead of $\boldsymbol{e_j}^*$.
