We study the spectral properties of a class of random matrices of the form $S_n^{-} = n^{-1}(X_1 X_2^* - X_2 X_1^*)$ where $X_k = Σ^{1/2}Z_k$, for $k=1,2$, $Z_k$'s are independent $p\times n$ complex-valued random matrices, and $Σ$ is a $p\times p$ positive semi-definite matrix, independent of the $Z_k$'s. We assume that $Z_k$'s have independent entries with zero mean and unit variance. The skew-symmetric/skew-Hermitian matrix $S_n^{-}$ will be referred to as a random commutator matrix associated with the samples $X_1$ and $X_2$. We show that, when the dimension $p$ and sample size $n$ increase simultaneously, so that $p/n \to c \in (0,\infty)$, there exists a limiting spectral distribution (LSD) for $S_n^{-}$, supported on the imaginary axis, under the assumptions that the spectral distribution of $Σ$ converges weakly and the entries of $Z_k$'s have moments of sufficiently high order. This nonrandom LSD can be described through its Stieltjes transform, which satisfies a coupled Marčenko-Pastur-type functional equations. In the special case when $Σ= I_p$, we show that the LSD of $S_n^{-}$ is a mixture of a degenerate distribution at zero (with positive mass if $c > 2$), and a continuous distribution with a symmetric density function supported on a compact interval on the imaginary axis. Moreover, we show that the companion matrix $S_n^{+} = Σ_n^\frac{1}{2}(Z_1Z_2^* + Z_2Z_1^*)Σ_n^\frac{1}{2}$, under identical assumptions, has an LSD supported on the real line, which can be similarly characterized.