A fully coupled perturbed formalism is developed fdr spin-coupled wave functions for the calculation of electric and magnetic response properties of molecules, particularly as functions of the nuclear geometry. The theory is specialized to the case of polarizabilities and magnetizabilities, but is easily generalized to include gradients and second derivatives with respect to nuclear coordinates and to spin-dependent properties. For polarizabilities and similar properties, for which the perturbation is purely real, the main quantity required is the Hessian matrix. This is in any case calculated in the stabilized Newton-Raphson method during the process of minimizing the energy. For magnetic properties, or for those properties for which the perturbation is purely imaginary, a modified Hessian is needed. This last is quickly assembled from the same quantities (first-, second-, third-, and fourth-order density matrices) that are computed for the normal Hessian.