The theory of relaxed density matrices has been developed for the calculation of second-order response properties at third-order Moller-Plesset (MP3) and linearized coupled cluster double excitation (L-CCD) levels of theory. The ensuing algorithm is applied to the determination of the isotropic and anisotropic nuclear magnetic shielding constants and magnetizabilities for thirteen molecules (H-2, N-2, F-2, HF, CO, KCN, HNC, H2O, NH3, H2O(2), HCHO, CH4, and HCCH). The method uses conventional gauge-variant orbitals but, by using large basis sets, produces results which are equivalent to those found with gauge-including orbitals. In general the L-CCD values of the manetizabilities are consistent with those obtained with multiconfigurational self-consistent-field (MCSCF) methods for molecules which have been treated by this method. For the nuclear shieldings, when there is a difference between L-CCD and MP3, the former gives results closer to the coupled-cluster singles and doubles level treatment augmented with a perturbation correction for connected triple excitations [CCSD(T)] which is our reference point. Again the results for the shieldings at the L-CCD level are quite good. We also use the paramagnetic components of the shieldings and magnetizabilities to determine the spin-rotation constants and rotational g tensors, respectively. These quantities are important since they may be compared more directly with experiment than the magnetizabilities and shieldings.