This paper numerically investigates particle saltation in a turbulent channel flow having a rough bed consisting of two to three layers of densely packed spheres. The Shields function is 0.065 which is just above the sediment entrainment threshold to give a bed-load regime. The applied methodology is a combination of three technologies, i.e., the direct numerical simulation of turbulent flow; the combined finite-discrete element modeling of the deformation, movement, and collision of the particles; and the immersed boundary method for the fluid-solid interaction. It is shown that the presence of entrained particles significantly modifies the flow profiles of velocity, turbulent intensities, and shear stresses in the vicinity of a rough bed. The quasi-streamwise-aligned streaky structures are not observed in the near-wall region and the particles scatter on the rough bed owing to their large size. However, in the outer flow region, the turbulent coherent structures recover due to the weakening rough-bed effects and particle interferences. First- and second-order statistical features of particle translational and angular velocities, together with sediment concentration and volumetric flux density profiles, are presented. Several key parameters of the particle saltation trajectory are calculated and agree closely with published experimental data. Time histories of the hydrodynamic forces exerted upon a typical saltating particle, together with those of the particle's coordinates and velocities, are presented. A strong correlation is shown between the abruptly decreasing streamwise velocity and increasing vertical velocity at collision which indicates that the continuous saltation of large-grain-size particles is controlled by collision parameters such as particle incident angle, local bed packing arrangement, and particle density, etc.