At least when $c:=x_i-y_i$ does not depend on $i$, any $n\times n$ submatrix of $A$ is of the form $M=(f_c(u_i-u_j)\colon i,j=1,\dots,n\}$, where $f_c(x):=f(x+c)$. So, $\det M>0$, since the function $f_c$ (being the characteristic function of a normal distribution) is positive definite.
However, $f_c=f_{c;2}$ is only "borderline" positive definite, in the sense that for any real $p>0$, the function $f_{c;p}$ defined by the formula $f_{c;p}(x):=\exp(-|(x+c)/\sigma|^p/2)$ is positive definite if and only if $p\le2$. So then, $\det M$ may be very close to $0$ and hence $A$ may be very close to a non-full rank matrix.
The same will hold if for some permutation $\pi$ of the set $[n]:=\{1,\dots,n\}$ we have that $x_i-y_{\pi(i)}$ does not depend on $i$. A similar effect can be expected if $x_i-y_{\pi(i)}$ does not much depend on $i$ for for some permutation $\pi$ of $[n]$. Because there a great number ($n!$) of such permutations, it is likely that there is one with this "does not much depend on $i$" property.
There also is the general tendency for large $n\times n$ matrices to be close to singular. Indeed, the determinant is the product of the eigenvalues (repeated according to their multiplicities). So, if even one eigenvalue is very close to $0$, then the determinant may be close to $0$. This effect is exponentially strong in the number of eigenvalues close to $0$.
That is why different regularization methods, such as Tikhonov regularization and Lasso method, are used with the method of least squares.