In this paper, we propose a robust minimum maximum mean square error Tomlinson-Harashima precoding (Min-Max-MSE THP) scheme and a low-complexity robust Min-Max-MSE hierarchical linear/THP (HL-THP) scheme for downlink massive multiuser multiple-input-multiple-output (MU-MIMO) systems with imperfect channel state information (CSI) at the transmitter. The proposed robust Min-Max-MSE HL-THP scheme comprises an inner linear beamformer (BF), which is designed based on second-order CSI statistics, and outer THP modules, which exploit the instantaneous overall CSI of the cascade of the actual channel and the inner BF. Thereby, the user terminals are divided into groups, where for each group a THP module successively mitigates the intra-group interference, whereas the inter-group interference is canceled by the inner BF. To ensure fairness, we adopt the maximization of the asymptotic signal-to-leakage-plus-noise ratio in the large system limit and the Min-Max-MSE as an optimization criterion for designing the inner BF and the per-group THP modules, respectively. Our analytical and simulation results show that the proposed robust Min-Max-MSE HL-THP scheme achieves a substantially improved performance in terms of the Max-MSE, maximum bit error rate, and minimum rate compared to linear regularized zero-forcing precoding. Moreover, the performance loss of the proposed robust Min-Max-MSE HL-THP scheme compared to the robust Min-Max-MSE THP scheme is small. In addition, our complexity analysis reveals that the proposed robust Min-Max-MSE HL-THP scheme has a much lower computational complexity than the Min-Max-MSE THP scheme. Hence, the robust Min-Max-MSE HL-THP scheme provides a favorable tradeoff between complexity and performance.