Different ways to incorporate two-dimensional systems which are not amenable to separation of variables treatment into the framework of supersymmetrical quantum mechanics (SUSY QM) are analysed. In particular, the direct generalization of one-dimensional Witten SUSY QM is based on supercharges of first order in momenta and allows one to connect the eigenvalues and eigenfunctions of two scalar and one matrix Schrödinger operators. The use of second-order supercharges leads to polynomial supersymmetry and relates a pair of scalar Hamiltonians, giving a set of partner systems with almost coinciding spectra. This class of systems can be studied by means of a new method of SUSY separation of variables, where supercharges allow separation of variables, but Hamiltonians do not. The method of shape invariance is generalized to two-dimensional models in order to construct purely algebraically a chain of eigenstates and eigenvalues for generalized Morse potential models in two dimensions.