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][1] and [Lasso method][2], are used with the method of least squares. 


  [1]: https://en.wikipedia.org/wiki/Least_squares#Tikhonov_regularization
  [2]: https://en.wikipedia.org/wiki/Least_squares#Lasso_method