To construct a N-body anti-symmetric wave function some derivations start with the requirement that the N-body wave function should be anti-symmetric under a permutation of coordinates, other derivations start with the requirement that the total wave function should be anti-symmetric under a permutation of labels or states, for example this derivation . I can see the equivalence is trivial if you assume you can freely change the order of the wavefunctions. This is true in the case the wavefunctions are just scalars but what if you are dealing with multicomponent wavefunctions. For example relativistic wavefunctions which are 4x1 matrices.
A simple example of my problem is the following, suppose we want to write down a general 2 body wavefunction in terms of single particle wavefunctions, \begin{align} \Psi(x_1,x_2) \propto \varphi_{a}(x_1)\varphi_b(x_2) \end{align} Where $x_i$ denote the spatial,spin,... coordinates, and $a,b$ denote single particle eigenstates. If we permute the coordinates in order to derive a antisymmetric wavefunction we get, \begin{align} \frac{1}{\sqrt{2}} \left(\varphi_{a}(x_1)\varphi_b(x_2) - \varphi_{a}(x_2)\varphi_b(x_1) \right) \end{align} while a permutation in states gives, \begin{align} \frac{1}{\sqrt{2}} \left(\varphi_{a}(x_1)\varphi_b(x_2) - \varphi_{b}(x_1)\varphi_a(x_2) \right) \end{align} Obviously above expressions are identical for scalar wavefunctions. When I try to check the normalisation of the two particle wavefunction (assuming the single particle states are orthormal and general $N \times 1$ matrices) I get for the label permutation, \begin{align} \int | \Psi(x_1,x_2) |^{2} &= \frac{1}{2} \int \varphi_{b}^{\dagger}(x_2) \varphi_{a}^{\dagger}(x_1) \varphi_{a}(x_1)\varphi_b(x_2) - \frac{1}{2} \int \varphi_{b}^{\dagger}(x_2) \varphi_{a}^{\dagger}(x_1)\varphi_{b}(x_1)\varphi_a(x_2) \\ &- \frac{1}{2} \int \varphi_{a}^{\dagger}(x_2) \varphi_{b}^{\dagger}(x_1) \varphi_{a}(x_1)\varphi_b(x_2) + \frac{1}{2} \int \varphi_{a}^{\dagger}(x_2) \varphi_{b}^{\dagger}(x_1)\varphi_{b}(x_1)\varphi_a(x_2) \\ &= \frac{1}{2} - 0 - 0 + \frac{1}{2} = 1 \end{align} Which is the expected result (note that the integral is assumed to integrate over all continuous degrees of freedom, or sum over the discrete ones). For the coordinate permutation however I get, \begin{align} \int | \Psi(x_1,x_2) |^{2} &= \frac{1}{2} \int \varphi_{b}^{\dagger}(x_2) \varphi_{a}^{\dagger}(x_1) \varphi_{a}(x_1)\varphi_b(x_2) - \frac{1}{2}\int \varphi_{b}^{\dagger}(x_2) \varphi_{a}^{\dagger}(x_1)\varphi_{a}(x_2)\varphi_b(x_1) \\ &- \frac{1}{2}\int \varphi_{b}^{\dagger}(x_1) \varphi_{a}^{\dagger}(x_2) \varphi_{a}(x_1)\varphi_b(x_2) + \frac{1}{2}\int \varphi_{b}^{\dagger}(x_1) \varphi_{a}^{\dagger}(x_2)\varphi_{a}(x_2)\varphi_b(x_1) \\ \end{align} Where it is not clear to my why the cross terms, \begin{align} \propto \int \textrm{d}x_1 \textrm{d}x_2 \left[ \varphi_{b}^{\dagger}(x_1) \varphi_b(x_2) \right] \left[ \varphi_{a}^{\dagger}(x_2) \varphi_{a}(x_1) \right] \end{align} should vanish. (Even if they did my main question still stands namely if you can always freely change the order of the wavefunctions). Note that the choice of permutation in coordinates or labels is equivalent to the choice whether to expand to rows or columns in a Slater determinant.
EDIT: upon request a more elaborate version. Suppose we want to write down a general wave function for 2 identical fermions $\Psi$. This wavefunction has to be antisymmetric. Suppose the single particle states are given by $$ \varphi_{\alpha_i}(r_i,s_i) \,\,\,\,\,\ i \in \{1,2\} $$ Where $\alpha_i$ denotes the quantum numbers of a state, while $r_i,s_i$ denote the spatial and spin coordinates. A way to antisymmeterize $N$-body wavefunction in terms of single particle wavefunctions is often introduced by the concept of a Slater determinant. In the case of two particles a Slater determinant looks like \begin{align} \Psi_{\alpha_1,\alpha_2}(r_1,s_1,r_2,s_2) = \frac{1}{\sqrt{2}} \left| \begin{array}{c c} \varphi_{\alpha_1}(r_1,s_1) & \varphi_{\alpha_2}(r_1,s_1) \\ \varphi_{\alpha_1}(r_2,s_2) & \varphi_{\alpha_2}(r_2,s_2) \end{array} \right|. \end{align} Now we can calculate this determinant by expanding along a row or along a column. Which is equivalent to the choice whether to permute the quantum numbers of the single particle states or the the coordinates of the single particle states. More explicitly, expanding the matrix along the first row (or equivalently, permuting the label of the single particle states) we get, $$ \frac{1}{\sqrt{2}} \left(\varphi_{\alpha_1}(r_1,s_1)\varphi_{\alpha_2}(r_2,s_2) - \varphi_{\alpha_2}(r_1,s_1)\varphi_{\alpha_1}(r_2,s_2) \right). $$ While an expansion along the first column gives (equivalent to a coordinate permutation), $$ \frac{1}{\sqrt{2}} \left(\varphi_{\alpha_1}(r_1,s_1)\varphi_{\alpha_2}(r_2,s_2) - \varphi_{\alpha_1}(r_2,s_2)\varphi_{\alpha_2}(r_1,s_1) \right). $$ The two above expression certainly look to be equivalent. However if we consider the wave functions to be $N$-component vectors I get some strange results. Suppose we want to calculate the normalization of the two particle wave function. We would then write the following $$ \sum_{s_1} \sum_{s_2} \int \textrm{d}r_1 \textrm{d}r_2 \Psi^{\dagger}(r_1,s_1,r_2,s_2) \Psi(r_1,s_1,r_2,s_2) $$ Note that I have introduced the $\dagger$ instead of the complex conjugate $*$ as the normalization should give a scalar. If we now replace $\Psi$ with the expression for the "label" permutation case we get, \begin{align*} \sum_{s_1} \sum_{s_2} & \int \textrm{d}r_1 \textrm{d}r_2 \Psi^{\dagger}(r_1,s_1,r_2,s_2) \Psi(r_1,s_1,r_2,s_2) \\ &= \frac{1}{2} \sum_{s_1} \sum_{s_2} \int \textrm{d}r_1 \textrm{d}r_2 \left(\varphi_{\alpha_1}(r_1,s_1)\varphi_{\alpha_2}(r_2,s_2) - \varphi_{\alpha_2}(r_1,s_1)\varphi_{\alpha_1}(r_2,s_2) \right)^{\dagger} \\ & \times \left(\varphi_{\alpha_1}(r_1,s_1)\varphi_{\alpha_2}(r_2,s_2) - \varphi_{\alpha_2}(r_1,s_1)\varphi_{\alpha_1}(r_2,s_2) \right) \\ &= \frac{1}{2} \sum_{s_1} \sum_{s_2} \int \textrm{d}r_1 \textrm{d}r_2 \left( \varphi^{\dagger}_{\alpha_2}(r_2,s_2) \varphi^{\dagger}_{\alpha_1}(r_1,s_1) \varphi_{\alpha_1}(r_1,s_1)\varphi_{\alpha_2}(r_2,s_2) \right. \\ & \,\, - \varphi^{\dagger}_{\alpha_2}(r_2,s_2) \varphi^{\dagger}_{\alpha_1}(r_1,s_1) \varphi_{\alpha_2}(r_1,s_1)\varphi_{\alpha_1}(r_2,s_2) \\ & \,\, - \varphi^{\dagger}_{\alpha_1}(r_2,s_2) \varphi^{\dagger}_{\alpha_2}(r_1,s_1) \varphi_{\alpha_1}(r_1,s_1)\varphi_{\alpha_2}(r_2,s_2) \\ & \,\, \left. + \varphi^{\dagger}_{\alpha_1}(r_2,s_2) \varphi^{\dagger}_{\alpha_2}(r_1,s_1) \varphi_{\alpha_2}(r_1,s_1)\varphi_{\alpha_1}(r_2,s_2) \right) \\ &= \frac{1}{2} \sum_{s_1} \int \textrm{d}r_1 \varphi^{\dagger}_{\alpha_1}(r_1,s_1) \varphi_{\alpha_1}(r_1,s_1) \sum_{s_2} \int \textrm{d}r_2 \varphi^{\dagger}_{\alpha_2}(r_2,s_2) \varphi_{\alpha_2}(r_2,s_2) \\ &-\frac{1}{2}\sum_{s_1} \int \textrm{d}r_1 \varphi^{\dagger}_{\alpha_1}(r_1,s_1) \varphi_{\alpha_2}(r_1,s_1) \sum_{s_2} \int \textrm{d}r_2 \varphi^{\dagger}_{\alpha_2}(r_2,s_2) \varphi_{\alpha_1}(r_2,s_2) \\ &-\frac{1}{2}\sum_{s_1} \int \textrm{d}r_1 \varphi^{\dagger}_{\alpha_2}(r_1,s_1) \varphi_{\alpha_1}(r_1,s_1) \sum_{s_2} \int \textrm{d}r_2 \varphi^{\dagger}_{\alpha_1}(r_2,s_2) \varphi_{\alpha_2}(r_2,s_2) \\ &+\frac{1}{2}\sum_{s_1} \int \textrm{d}r_1 \varphi^{\dagger}_{\alpha_2}(r_1,s_1) \varphi_{\alpha_2}(r_1,s_1) \sum_{s_2} \int \textrm{d}r_2 \varphi^{\dagger}_{\alpha_1}(r_2,s_2) \varphi_{\alpha_1}(r_2,s_2) \end{align*} Note that we have used the general expression $(AB)^{\dagger} = B^{\dagger} A^{\dagger}$ as well as the fact that the product $\varphi^{\dagger} \varphi$ contracts to a scalar and can be safely dragged trough the expressions. Now if one assumes the single particle states are orthogonal we have the following relation, \begin{align} \sum_{s_i}\int \textrm{d} r_i \varphi^{\dagger}_{\alpha_m}(r_i,s_i) \varphi_{\alpha_n}(r_i,s_i) = \delta_{m,n} \end{align} So we get, \begin{align} \sum_{s_1} \sum_{s_2} \int \textrm{d}r_1 \textrm{d}r_2 &\Psi^{\dagger}(r_1,s_1,r_2,s_2) \Psi(r_1,s_1,r_2,s_2) \\ & =\frac{1}{2} \delta_{\alpha_1,\alpha_1} \delta_{\alpha_2,\alpha_2} - \frac{1}{2} \delta_{\alpha_1,\alpha_2} \delta_{\alpha_2,\alpha_1} - \frac{1}{2} \delta_{\alpha_2,\alpha_1} \delta_{\alpha_1,\alpha_2} + \frac{1}{2} \delta_{\alpha_2,\alpha_2} \delta_{\alpha_1,\alpha_1} \\ &=1 \end{align} Which is to be expected. If we start from the second expression for $\Psi$ however we get the following, \begin{align*} \sum_{s_1} \sum_{s_2} & \int \textrm{d}r_1 \textrm{d}r_2 \Psi^{\dagger}(r_1,s_1,r_2,s_2) \Psi(r_1,s_1,r_2,s_2) \\ &= \frac{1}{2} \sum_{s_1} \sum_{s_2} \int \textrm{d}r_1 \textrm{d}r_2 \left(\varphi_{\alpha_1}(r_1,s_1)\varphi_{\alpha_2}(r_2,s_2) - \varphi_{\alpha_1}(r_2,s_2)\varphi_{\alpha_2}(r_1,s_1) \right)^{\dagger} \\ & \times \left(\varphi_{\alpha_1}(r_1,s_1)\varphi_{\alpha_2}(r_2,s_2) - \varphi_{\alpha_1}(r_2,s_2)\varphi_{\alpha_2}(r_1,s_1) \right) \\ &= \frac{1}{2} \sum_{s_1} \sum_{s_2} \int \textrm{d}r_1 \textrm{d}r_2 \left( \varphi^{\dagger}_{\alpha_2}(r_2,s_2) \varphi^{\dagger}_{\alpha_1}(r_1,s_1) \varphi_{\alpha_1}(r_1,s_1)\varphi_{\alpha_2}(r_2,s_2) \right. \\ & \,\, - \varphi^{\dagger}_{\alpha_2}(r_2,s_2) \varphi^{\dagger}_{\alpha_1}(r_1,s_1) \varphi_{\alpha_1}(r_2,s_2)\varphi_{\alpha_2}(r_1,s_1) \\ & \,\, - \varphi^{\dagger}_{\alpha_2}(r_1,s_1) \varphi^{\dagger}_{\alpha_1}(r_2,s_2) \varphi_{\alpha_1}(r_1,s_1)\varphi_{\alpha_2}(r_2,s_2) \\ & \,\, \left. + \varphi^{\dagger}_{\alpha_2}(r_1,s_1) \varphi^{\dagger}_{\alpha_1}(r_2,s_2) \varphi_{\alpha_1}(r_2,s_2)\varphi_{\alpha_2}(r_1,s_1) \right) \\ \end{align*} It is clear that the first and the fourth term will give $1/2 + 1/2 = 1$ as they are exactly the same as in the previous expression. However I do not see why the cross term should be necessarily zero. So it would seem that starting from the different antisymmeterization expression one gets different results. Hence my phrase "order matters" as the following appears to be true $$ \frac{1}{\sqrt{2}} \left(\varphi_{\alpha_1}(r_1,s_1)\varphi_{\alpha_2}(r_2,s_2) - \varphi_{\alpha_2}(r_1,s_1)\varphi_{\alpha_1}(r_2,s_2) \right) \neq \frac{1}{\sqrt{2}} \left(\varphi_{\alpha_1}(r_1,s_1)\varphi_{\alpha_2}(r_2,s_2) - \varphi_{\alpha_1}(r_2,s_2)\varphi_{\alpha_2}(r_1,s_1) \right) $$ A final note about the Dirac notation. If one would start from the Dirac formalism and introduce $ | \Psi \rangle$ as $$ | \Psi \rangle = \frac{1}{\sqrt{2}} \left( | \alpha_1 \rangle | \alpha_2 \rangle - | \alpha_2 \rangle | \alpha_1 \rangle \right) $$ A projection into coordinate space gives the "label" permutation expression, $$ \langle r_1,s_1 ; r_2,s_2 | \Psi \rangle = \frac{1}{\sqrt{2}} \left( \langle r_1,s_1| \alpha_1 \rangle \langle r_2,s_2 | \alpha_2\rangle - \langle r_1,s_1| \alpha_2 \rangle \langle r_2,s_2 | \alpha_1\rangle \right) \\ = \frac{1}{\sqrt{2}} \left(\varphi_{\alpha_1}(r_1,s_1)\varphi_{\alpha_2}(r_2,s_2) - \varphi_{\alpha_2}(r_1,s_1)\varphi_{\alpha_1}(r_2,s_2) \right) $$ The expression for the normalization starting from the Dirac formalism gives, $$ \langle \Psi| \Psi \rangle = \frac{1}{2} \left( \langle \alpha_1 | \alpha_1 \rangle \langle \alpha_2 | \alpha_2 \rangle - \langle \alpha_1 | \alpha_2 \rangle \langle \alpha_2 | \alpha_1 \rangle - \langle \alpha_2 | \alpha_1 \rangle \langle \alpha_1 | \alpha_2 \rangle + \langle \alpha_2 | \alpha_2 \rangle \langle \alpha_1 | \alpha_1 \rangle \right) $$ Using the unit identities: $1 = \sum_s | s \rangle \langle s |$ and, $1 = \int \textrm{d} r | r \rangle \langle r |$ we get, $$ \langle \Psi| \Psi \rangle = \frac{1}{2} \sum_{s_1} \sum_{s_2} \int \textrm{d} r_1 \int \textrm{d} r_2 \\ \left( \langle \alpha_1 | r_1,s_1 \rangle \langle r_1,s_1 | \alpha_1 \rangle \langle \alpha_2 | r_2,s_2 \rangle \langle r_2,s_2| \alpha_2 \rangle - \langle \alpha_1 | r_1,s_1 \rangle \langle r_1,s_1 | \alpha_2 \rangle \langle \alpha_2 | r_2,s_2 \rangle \langle r_2,s_2| \alpha_1 \rangle - \langle \alpha_2 | r_1,s_1 \rangle \langle r_1,s_1| \alpha_1 \rangle \langle \alpha_1 | r_2,s_2 \rangle \langle r_2,s_2| \alpha_2 \rangle + \langle \alpha_2 | r_1,s_1 \rangle \langle r_1,s_1| \alpha_2 \rangle \langle \alpha_1 | r_2,s_2 \rangle \langle r_2,s_2| \alpha_1 \rangle \right) $$ Which is when using $ \langle r_i, s_i | \alpha_i \rangle = \varphi_{\alpha_i}(r_i,s_i)$ exactly the same as the normalization expression I get for the label permutation case, which works out nicely to $1$ as expected.