In this paper, we focus on the numerical simulation for two-dimensional nonlinear fractional sub-diffusion equations in the presence of time delay. Firstly, we investigate the existence, uniqueness and regularity of the solution for such problems. The theoretical result implies that the solution at τ+ is smoother than that at 0+, where τ is a constant time delay, and this is an improvement for the work (Tan et al., 2022). Secondly, a high-order difference scheme based on L1 method is constructed. For the sake of repairing the convergence order in temporal direction and improving the computational efficiency, an efficient time two-grid algorithm based on nonuniform meshes is first developed. The convergence order of the two-grid scheme reaches O(NF−min{rα,2−α}+NC−min{2rα,4−2α}+h12+h22), where NF and NC represent the number of the fine and coarse grids respectively, while h1 and h2 are the space-step sizes. Furthermore, stability and convergence analysis of the proposed scheme are carefully verified by energy method. Finally, numerical experiments are carried out to show the validity of theoretical statements.