I'll give an alternative derivation. Let me just say that I will work entirely in Euclidean signature for simplicity. Generalizing to Lorentzian just requires an extra $\det(\eta_{\mu\nu}) = (-1)^{n-1}$ overall.
By definition of determinant, given an $n\times n$ matrix $A_{ij}$ one has
$$
\det A = \epsilon_{i_1\cdots i_n} A_{1,i_1}\cdots A_{n,i_n}\,.\label{1}\tag{1}
$$
Here we will use Einstein convention throughout. This is also equivalent to
$$
\det A = \frac{1}{n!}\epsilon_{i_1\cdots i_n}\epsilon_{j_1\cdots j_n} A_{j_1,i_1}\cdots A_{j_n,i_n}\,.\tag{2}\label{2}
$$
See the part at the botton if you want a proof of it. At the same time we have
$$
\begin{aligned}
\epsilon_{k_1\cdots k_n}\det A &= \epsilon_{k_1\cdots k_n}\epsilon_{i_1\cdots i_n} A_{1,i_1}\cdots A_{n,i_n} =\\&=
\epsilon_{i_1\cdots i_n} A_{k_1,i_1}\cdots A_{k_n,i_n}\,.
\end{aligned}\tag{3}\label{3}
$$
The way to show \eqref{3} is by sorting the $A_{a,i_a}$ so that they appear with the row indices matching $k_a$, then we apply the inverse permutation to the $i_a$ to put them back in the initial order (call this permutation $\sigma$). This comes at a cost of introducing a $\mathrm{sgn}(\sigma)$, i.e. the parity of the permutation, which is precisely $\epsilon_{k_1\cdots k_n}$, thus cancelling that factor (since either $1^2$ or $(-1)^2$ is $1$).
Just as a sanity check: contract \eqref{3} with $\epsilon_{k_1\cdots k_n}$ and compare with \eqref{2} to see that, indeed
$$
\epsilon_{k_1\cdots k_n}\epsilon_{k_1\cdots k_n} = n!\,.
$$
But back to business. I will now do something weird. Let me take the entries of $A$ to be symbols, namely
$$
A_{i,j} \equiv \delta_{i, K_j}\,.
$$
By that I do not mean the identity matrix but rather
$$
A = \left(\begin{matrix}
\delta_{1,K_1} &\delta_{1,K_2} &\ldots&\delta_{1,K_n}\\
\delta_{2,K_1} &\delta_{2,K_2} \\
\vdots && \ddots\\
\delta_{n,K_1}&&&\delta_{n,K_n}
\end{matrix}\right)\,.
$$
The symbols still commute, so everything goes through, but now \eqref{2} and \eqref{3} say
$$
\epsilon_{i_1\cdots i_n} \delta_{k_1,K_{i_1}}\cdots \delta_{k_n,K_{i_n}} = \frac1{n!} \epsilon_{k_1\cdots k_n}\epsilon_{i_1\cdots i_n}\epsilon_{j_1\cdots j_n}\delta_{j_1,K_{i_1}}\cdots \delta_{j_n,K_{i_n}}\,.
$$
That's a lot of indices, I apologize. But we can work our way to a more readable expression. First let's just contract the $\delta$'s in the right hand side
$$
\epsilon_{i_1\cdots i_n} \delta_{k_1,K_{i_1}}\cdots \delta_{k_n,K_{i_n}} = \frac1{n!} \epsilon_{k_1\cdots k_n}\epsilon_{i_1\cdots i_n}\epsilon_{K_{i_1}\cdots K_{i_n}}\,.
$$
Now notice that the sum over the $i_a$'s in the RHS produces $n!$ times the same term, because I can either have an even or odd permutation, but the sign of the $\epsilon_{i_1\ldots}$ tensor is the same as that of the $\epsilon_{K_{i_1}\cdots}$ tensor. So
$$
\epsilon_{i_1\cdots i_n} \delta_{k_1,K_{i_1}}\cdots \delta_{k_n,K_{i_n}} = \epsilon_{k_1\cdots k_n}\epsilon_{K_1\cdots K_n}\,.
$$
Finally, the left hand side is by definition the antisymmetrization of the $\delta$'s over the second index. By convention the antisymmetrization has weight $1$ (see here or here), meaning that it doesn't overcount. Since here we instead have $n!$ terms we have to multiply by that.
$$
n!\,\delta_{k_1,[K_1}\cdots \delta_{k_n,K_n]} = \epsilon_{k_1\cdots k_n}\epsilon_{K_1\cdots K_n}\,.
$$
I also replaced $i_a$ by $a$ because they are antisymmetrized so the order they had initially doesn't matter and I can re-sort them as I please (paying signs of course).
Proof of the equality \eqref{1} \eqref{2}.
Take the expression \eqref{1} and exchange the position of $A_{p, i_p}$ and $A_{q, i_q}$. For the sake of concreteness let's say $p<q$. This does nothing because the entries of the matrix are just numbers!
$$
\begin{aligned}
\det A &= \epsilon_{i_1\cdots i_n} A_{1,i_1}\cdots A_{p, i_p}\cdots A_{q, i_q}\cdots A_{n,i_n}
\\& = \epsilon_{i_1\cdots i_n} A_{1,i_1}\cdots A_{q, i_q}\cdots A_{p, i_p}\cdots A_{n,i_n}\,.
\end{aligned}
$$
Ok, nothing happened, but let me now swap $i_p$ and $i_q$ in the $\epsilon$:
$$
\epsilon_{i_1\cdots i_p\cdots i_q\cdots i_n} = -\epsilon_{i_1\cdots i_q\cdots i_p\cdots i_n}\,.
$$
Obvious! We get a minus sign. So therefore I'll just rename $i_p$ to $i_q$ and $i_q$ to $i_p$ (I'm always free to do so since there are summed over)
\begin{aligned}
\det A &= \epsilon_{i_1\cdots i_n} A_{1,i_1}\cdots A_{p, i_p}\cdots A_{q, i_q}\cdots A_{n,i_n}
\\& = -\epsilon_{i_1\cdots i_n} A_{1,i_1}\cdots A_{q, i_p}\cdots A_{p, i_q}\cdots A_{n,i_n}\,.
\end{aligned}
We have just proven that in the product $A_{1,i_1}\cdots A_{n,i_n}$ we can antisymmetrize over the row indices as well. That's because, as we saw, swapping any two row indices gives the same contribution up to a sign. We can then say
$$
\det A = \frac{1}{n!}\sum_{\sigma \in S_n}\mathrm{sgn}(\sigma)\,\epsilon_{i_1\cdots i_n} A_{\sigma(1),i_1}\cdots A_{\sigma(n),i_n}\,,
$$
where $S_n$ is the permutation group of $n$ elements and $\mathrm{sgn}(\sigma)$ is the parity of the permutation. I divided by $n!$ because every term is now counted $|S_n|= n!$ times. You might know that, in general, for any tensor
$$
\sum_{\sigma \in S_n}\mathrm{sgn}(\sigma) \,T_{\sigma(1)\cdots \sigma(n)} = \epsilon_{i_1\ldots i_n} T_{i_1\ldots i_n}\,.
$$
This is basically by definition of $\epsilon_{i_1\ldots i_n}$. Looking back, we just proved \eqref{2}.