The self-consistent random phase approximation (RPA), based on the framework of relativistic energy density functionals, is employed in the study of isovector and isoscalar dipole response in $^{68}$Ni, $^{132}$Sn, and $^{208}$Pb. The evolution of pygmy dipole states (PDS) in the region of low excitation energies is analyzed as a function of the density-dependence of the symmetry energy for a set of relativistic effective interactions. The occurrence of PDS is predicted in the response to both the isovector and isoscalar dipole operators, and its strength is enhanced with the increase of the symmetry energy at saturation and the slope of the symmetry energy. In both channels the PDS exhausts a relatively small fraction of the energy-weighted sum rule but a much larger percentage of the inverse energy-weighted sum rule. For the isovector dipole operator the reduced transition probability $B(E1)$ of the PDS is generally small because of pronounced cancellation of neutron and proton partial contributions. The isoscalar reduced transition amplitude is predominantly determined by neutron particle-hole configurations, most of which add coherently, and this results in a collective response of the PDS to the isoscalar dipole operator.