We present the results of the first quantum chemical investigations of H-1 NMR hyperfine shifts in the blue copper proteins (BCPs): amicyanin, azurin, pseudoazurin, plastocyanin, stellacyanin, and rusticyanin. We find that very large structural models that incorporate extensive hydrogen bond networks, as well as geometry optimization, are required to reproduce the experimental NMR hyperfine shift results, the best theory vs experiment predictions having R-2 = 0.94, a slope = 1.01, and a SD = 40.5 ppm (or similar to 4.7% of the overall similar to 860 ppm shift range). We also find interesting correlations between the hyperfine shifts and the bond and ring critical point properties computed using atoms-in-molecules theory, in addition to finding that hyperfine shifts can be well-predicted by using an empirical model, based on the geometry-optimized structures, which in the future should be of use in structure refinement.