The direct (recomputation of two-electron integrals) implementation of the gauge-including atomic orbital (GIAO) and the CSGT (continuous set of gauge transformations) methods for calculating nuclear magnetic shielding tensors at both the Hartree-Fock and density functional levels of theory are presented. Isotropic C-13, N-15, and O-17 magnetic shielding constants for several molecules, including taxol (C47H51NO14 using 1032 basis functions) are reported. Shielding tensor components determined using the GIAO and CSGT methods are found to converge to the same value at sufficiently large basis sets; however, GIAO shielding tensor components for atoms other than carbon are found to converge faster with respect to basis set size than those determined using the CSGT method for both Hartree-Fock and DFT. For molecules where electron correlation effects are significant, shielding constants determined using (gradient-corrected) pure DFT or hybrid methods (including a mixture of Hartree-Fock exchange and DFT exchange-correlation) are closer to experiment than those determined at the Hartree-Fock level of theory. For the series of molecules studied here, the RMS error for C-13 chemical shifts relative to TMS determined using the B3LYP hybrid functional with the 6 - 311 + G(2d,p) basis is nearly three times smaller than the RMS error for shifts determined using Hartree-Fock at this same basis. Hartree-Fock C-13 chemical shifts calculated using the 6 - 31G* basis set give nearly the same RMS error as compared to experiment as chemical shifts obtained using Hartree-Fock with the bigger 6 - 311 + G(2d,p) basis set for the range of molecules studied here. The RMS error for chemical shifts relative to TMS calculated at the Hartree-Fock 6 - 31G* level of theory for taxol (C47H51NO14) is 6.4 ppm, indicating that for large systems, this level of theory is sufficient to determine accurate C-13 chemical shifts.