Fifteen procedures based on hybrid density functional theory were used to calculate magnetic properties for the carbon-bound hydrogen nuclei of SO small to modest-sized organic molecules. The predicted isotropic shieldings derived from the various methods were compared with each other and also with solution experimental data. The computational methods investigated included the IGAIM and GIAO procedures, the 6-311++G(d,p), 6-311++G(2df,p), and 6-311++G(3df,2p) basis sets, the B3LYP, B3P86, and B3PW91 hybrid density functionals, and molecular geometries optimized using both MP2 and B3LYP methods. Although agreement with experiment consistently improved as the basis set was enlarged, the improvement upon going from 6-311++G(2df,p) to 6-311++G(3df,2p) was extremely small, and even the difference between 6-311++G(d,p) and 6-311++G(2df,p) was of a. modest size. The CIAO and IGAIM procedures yielded very similar results in conjunction with the largest basis set, but CIAO suffered considerably less degradation than did IGAIM as the basis set size was decreased. The three functionals B3LYP, B3P86, and B3PW91 performed in al extremely similar fashion, although B3LYP proved marginally superior to the others. The method of geometry optimization also was found to make little difference. Of the computational methods investigated, the GIAO/B3LYP/6-311++G(d,p)//B3LYP/6-31+G(d) procedure probably represents the best compromise between accuracy and expense and yielded proton chemical shifts having a root-mean-square error of 0.15 ppm in comparison with solution experimental values after empirical linear scaling. The more expensive GIAO/B3LYP/6-311++G(2df,p)//B3LYP/6-31+G(d) method provided only a slightly lower root-mean-square error of 0.14 ppm.